LHC lifetime frontier and visible decay searches in composite asymmetric dark matter models

The LHC lifetime frontier will probe dark sector in near future, and the visible decay searches at fixed-target experiments have been exploring dark sector. Composite asymmetric dark matter with dark photon portal is a promising framework explaining the coincidence problem between dark matter and visible matter. Dark strong dynamics provides rich structure in the dark sector: the lightest dark nucleon is the dark matter, while strong annihilation into dark pions depletes the symmetric components of the dark matter. Dark photons alleviate cosmological problems. Meanwhile, dark photons make dark hadrons long-lived in terrestrial experiments. Moreover, the dark hadrons are produced through the very same dark photon. In this study, we discuss the visible decay searches for composite asymmetric dark matter models. For a few GeV dark nucleons, the LHC lifetime frontier, MATHUSLA and FASER, has a potential to discover their decay when kinetic mixing angle of dark photon is ϵ ≳ 10−4. On the other hand, fixed-target experiments, in particular SeaQuest, will have a great sensitivity to dark pions with a mass below GeV and with kinetic mixing ϵ ≳ 10−4 in addition to the LHC lifetime frontier. These projected sensitivities to dark hadrons in dark photon parameter space are comparable with the future sensitivities of dark photon searches, such as Belle-II and LHCb.


Introduction
Lifetime is one of the most important properties of particles. In the Standard Model (SM), some particles are long-lived: for example, proton has a very long lifetime of 10 34 years, neutron decays into protons with a lifetime of O(10 3 ) s, and charged pions has a lifetime of O(10) ns. As in the SM sector, a dark sector where dark matter (DM) particles reside has possibly a rich structure: it is very natural that dark particles have a variety of lifetime, and some dark particles may be long-lived. At least, DM particles should have a lifetime longer than the age of the Universe: this is one of the known properties of DM particles.
The last decade has seen efforts in searching for the decay signal of long-lived particles in a dark sector by use of data collected at the past fixed-target experiments [1][2][3][4][5][6][7][8][9]. The main target particles in these experiments should have a lifetime of O(1) µs. The darksector particles may scatter with the SM particles via mediator particles that connect the dark sector with the SM sector. The scattering signals have been also explored at fixedtarget experiments with downstream detectors such as LSND [10,11], E137 [12,13], and MiniBooNE [14]. The produced dark-sector particles may decay into lighter dark states, which is referred to as invisible decay search: for instance, the mono-photon search at BaBar [15] puts a constraint on invisible decay of dark photon. Besides reanalyzing data from the existing experiments, dark sector models will be tested at the lifetime frontier by currently running experiments, upgrades of the existing experiments, and newly proposed experiments: displaced vertex searches at HPS [16], SHiP [17,18], and SeaQuest [19][20][21], and invisible decay at LDMX [22,23] and Belle-II [24,25]. These studies have been limited to the searches for the long-lived particle with mass of GeV or below since fixed-target experiments utilize low-energy beams.
The LHC experiment utilizes more energetic beams, and hence it is possible to produce long-lived particles with the mass of GeV or above. Various strategies have been proposed to search for the long-lived particles at the existing LHC experiments near ATLAS, CMS, and LHCb detectors [26][27][28]: they are the displaced vertex searches, and the future sensitivity of the long-lived particles is limited by the detector sizes and the tracking system of the LHC experiments. In particular, since there is no shield to reduce background events in the experiments, it is important to reconstruct the decay vertices in the detector. Therefore, lifetime of interest is O(10-10 2 ) ps in the LHC experiments with precise vertex reconstruction. Meanwhile, a variety of experiments at the LHC has been proposed to search for long-lived particles recently by locating additional detectors far from the interaction points of the LHC: MATHUSLA [29], FASER [30], and CODEX-b [31], collectively called LHC lifetime frontier. The detector locations are O(10-10 2 ) m away from the interaction point, and thus they would have a potential to explore the long-lived particles with a lifetime of O(10-10 2 ) ns.
We suggest that the composite asymmetric dark matter (ADM) framework, where DM particle masses are to be GeV and above, are very good target of the LHC lifetime frontier. The basic concept of the composite ADM framework is the following. Once the particle-antiparticle asymmetry is generated in either the SM or the dark sectors, it is shared among the two sectors via some portal interactions. Consequently asymmetries of -2 -

JHEP03(2022)176
baryonic matter and DM have the same origin and close to each other. Meanwhile, the observed energy densities of baryonic matter and DM are close to each other, Ω DM 5Ω B . Thus, DM particle masses are to be GeV. Compositeness plays three key roles: generation of the GeV mass, stability of the DM, and depletion of the symmetric DM component  (see ref. [78] for a review). The strong dynamics in the dark sector dynamically explains the DM mass in the GeV range through the dimensional transmutation. The accidental dark baryon number conservation ensures the longevity of the lightest dark nucleon, i.e., DM. Annihilation of the dark nucleon and anti-nucleon into the dark pions leaves only the asymmetric component, so that the DM asymmetry determines the DM relic density today as in the SM nucleons.
In this article, we specifically consider a composite ADM model with dark photon [67]. Dark photon plays a significant role in the composite ADM framework. Without dark photon, the dark sector contains a large entropy density since in the early Universe two sectors are in thermal contact via portal interaction sharing the generated asymmetry between two sectors, and cause cosmological problems. The dark photon decays into e + e − via the kinetic mixing with the SM photon [79], and it releases the entropy of the dark sector into the SM sector. We consider dark photons are the lightest particle in the dark sector in this study, and therefore the kinetic mixing of interest is below 10 −3 , evading constraints from searches for the prompt decay of dark photons at electron-positron collider experiments [3,8,9,[80][81][82] (see also review [83]).
Through the very same dark photon, dark protons that are charged under dark quantum electrodynamics (QED) are produced in collider experiments. Besides, through the dark photon, dark protons decay into dark neutrons that are neutral dark nucleons plus SM particles when dark protons are heavier than dark neutrons, and this transition can have a long lifetime due to the smallness of the kinetic mixing. The lightest dark nucleon must have the multi-GeV mass to be the ADM, and hence it is worthwhile for dark nucleon searches to explore the LHC lifetime frontier. In addition to the dark nucleon searches at the LHC frontier, a rich structure of the dark sector allows us to take advantage of various search strategies for the dark sector. When dark pions are heavier than dark photon but lighter than twice of the dark photon mass, dark pions decay into a dark photon emitting the SM particles via the anomaly-induced interaction. Dark pions are lighter than dark nucleons, and thus fixed-target experiments are also suitable for the dark pion visible decay signals in addition to the LHC lifetime frontier. At the same time, as we mentioned, dark photon has been searched for by visible decay searches at fixed-target experiments and prompt decay searches at collider experiments.
In this study, we investigate the sensitivity to the composite ADM models at the visible decay searches near future. We have found that the visible decay searches for the dark hadrons will explore the dark photon parameter space that the prompt decay searches of dark photon will also investigate. In particular, the LHC lifetime frontier will have a sensitivity to the dark nucleons with multi-GeV mass, and the sensitivity area corresponds to the dark photon parameter with a few GeV mass and with the kinetic mixing 10 −4 . As for dark pions, they are produced at both LHC and the fixed-target experiments. SeaQuest, one of the on-going fixed-target experiments, will particularly have a great sensitivity to the dark pion with sub-GeV and the kinetic mixing 10 −5 . The sensitivity area of the LHC lifetime frontier corresponds to the dark photon parameter with the sub-GeV mass and with the kinetic mixing 10 −4 . This paper is organized as follows. In section 2, we briefly review a composite model of ADM [67]. We discuss mass spectrum and lifetime of dark hadrons, and existing constraints except for the visible decay searches in section 3. In section 4, we discuss the visible decay searches of dark hadrons, and then we summarize sensitivities of visible decay searches in our composite ADM model in section 5. Section 6 is devoted to conclusions of our work.

Composite asymmetric dark matter model
We review the composite ADM scenario with a dark photon portal and fix our notation in this section. We consider the vector-like two-flavor SU(3) D × U(1) D dynamics proposed by ref. [67]. We list the minimal particle contents in the dark sector in table 1. The dark quarks are in (anti-)fundamental representation of SU(3) D , are charged under the dark QED U(1) D , and carry the same B − L number as the SM quarks. We refer to U and D as up and down dark quarks, respectively, since they have similar properties as the SM up and down quarks. The dark quarks are confined into dark hadrons below the dynamical scale Λ QCD of the dark QCD, SU(3) D . DM consists of the dark nucleons, whose stability is ensured by the B − L number conservation. The DM asymmetry has the same origin as the SM baryon asymmetry, and hence the mass of the dark nucleons is close to that of the SM nucleons as indicated from the coincidence of the mass densities of the baryonic matters and the DM, Ω DM ∼ 5Ω B .
Dark pions also play a significant role in the composite ADM models; the annihilation into the dark pions depletes the symmetric part of the DM. If the lightest pion is stable, it could lead to the overclosure of the Universe or a too large effective number of neutrino species. The dark QED, U(1) D , is therefore introduced to release the entropy of the dark sector to the SM sector. The dark photon gets its mass via spontaneous symmetry breaking of U(1) D by the dark Higgs φ D .
The Lagrangian density of the dark sector is given by

JHEP03(2022)176
Here, L gauge denotes the kinetic terms of the dark photon and the dark gluons, and L φ D contains the kinetic term and the potential term of the dark Higgs. The dark photon mass arises from the vacuum expectation value (VEV) v D of the dark Higgs, m 2 A = 2e 2 v 2 D with the U(1) D coupling e .
L portal connects the dark sector to the SM sector. This includes two kinds of portal interactions: one is called as the intermediate-scale portal interaction sharing the generated asymmetry between the SM sector and the dark sector, another is the dark photon portal interaction that releases the entropy in the dark sector to the SM sector. As for the intermediate-scale portal interaction, the charge assignment of the dark quarks in table 1 allows us to write the following interactions between our sector and the dark sector [67]: where F µν and F µν are field strength tensors of the SM photon and the dark photon, respectively. After diagonalizing the kinetic term, the dark photon A couples to the SM particles, and the coupling constant is proportional to their charge, eQ f , while the darksector particles remain electrically neutral. The dark photon decays into the SM particles through the coupling, where R( √ s) = σ(e + e − → hadrons)/σ(e + e − → µ + µ − ) takes into account the hadronic resonances [89,90]. The darkly charged particles interact with the SM particles through -5 -JHEP03(2022)176 dark photon exchange, and hence DM interact with ordinary matter when the DM is darkly charged.
The Lagrangian density L Q contains the dark quarks: Here, we use the hatted-notation for fields in the charge basis. The vector-like masses, m U and m D , are assumed to be smaller than the dynamical scale Λ QCD . The dark QCD interactions preserve the flavors of dark quarks, and hence there would be approximate global symmetry interchanging up and down dark quarks, which we refer to as isospin symmetry of the dark quarks. We introduce the isospin doublet of the dark quarks aŝ Q → e iα/2 e iβ/2 LQ ,Q → e iα/2 e −iβ/2Q R † , (2.6) where α and β are phases of U(1) A and U(1) B , while L and R denote SU(2) L and SU(2) R rotations, respectively. The chiral symmetry in the dark sector, SU(2) L × SU(2) R , is explicitly broken by the quark mass terms and U(1) D interactions. We systematically include the explicit breaking of chiral symmetry in the dark hadron spectrum discussed later by treating the breaking terms as spurion fields. The mass matrix of the dark quarks, M Q , has off-diagonal entries when the dark Higgs obtains the VEV v D : 1 The mass matrix transforms as a spurion under the global flavor symmetry as M Q → RM Q L † . This matrix is diagonalized by unitary matrices U L and U R as In the following, we assume that M Q is hermitian (parity-conserving), i.e., y 1 = y 2 , and thus the diagonalizing matrices are the same, U L = U R ≡ U V ; we also discuss the case that M Q is non-hermitian in appendix A. The covariant derivative and the charge matrix for dark quarks is defined by Here, we omit the dark QCD interaction. The charge matrices transform as spurions under the global flavor symmetry asq L → Lq L L † andq R → Rq R R † .

JHEP03(2022)176
The dark quarks in the charge basis and in the mass basis are related with each other by the unitary transformation: whereQ andQ indicate dark quarks in the charge basis, again, while Q and Q are the quarks in the mass basis. The charge matrices in the mass basis are given by 11) and are different as U L = U R . We discuss the long-lived dark hadrons in the composite ADM scenario in this article. The isospin symmetry is broken by the dark quark masses and the U(1) D interaction. This leads to the mass difference among states in an isospin multiplet even in hadronic picture below the dynamical scale, and thus the heavier states can decay into the lighter state. The heavier state can be a long-lived particle when the isospin symmetry is slightly broken.
We use the linear sigma model (LSM) as a low-energy effective theory to describe the mass spectrum of the dark hadrons and the interactions among the dark hadrons below Λ QCD . We discuss the LSM in detail in appendix A. The LSM field Φ is a matrix-form scalar field, and corresponds to the quark bilinear Φ ij ∼ Q j Q i where the dark quarks Q i and Q i are in mass basis. On the other hand, the isospin-doublet dark nucleons are denoted by the Weyl fermions N and N . The transformation of Φ , N , and N under the flavor symmetry is given by Here, L , R and α , β characterize the same transformation under the flavor symmetry as dark quarks. The low-energy effective Lagrangian consists of the LSM field Φ, the dark nucleons, and the dark photon. The LSM Lagrangian has the following form: Here, a VEV of Φ develops since µ 2 is positive, and hence the VEV breaks the global flavor symmetry SU(2) L × SU(2) R .
The covariant derivative of the LSM field is defined by where q L and q R are the charge matrices. We note that the charge matrices are identical (q L = q R ) when Φ is in the charge basis. The naïve dimensional analysis (NDA) with -7 -

JHEP03(2022)176
large-N C scaling of the low-energy parameters in the LSM is given by (see refs. [91,92] for NDA and refs. [93][94][95][96][97][98] for large-N C scaling) with the number of dark colors, N C = 3, and the chiral symmetry breaking scale in the dark sector, Λ χSB . The fourth term of V LSM comes from the instanton-induced quark interaction and breaks U(1) A [99,100], which we refer to as the anomaly term. We set the coefficient of the term c to be positive so that the term is parity symmetric, which directly gives the mass difference between the dark η meson and the dark pions in the two-flavor case. The last term originates from the quark mass term, and thus H is hermition, and H ≡ j a T a = j 0 T 0 +j 3 T 3 is proportional to the diagonalized mass matrix of the dark quarks, We define scalar and pseudoscalar components of Φ as follows.
where the matrixλ 0 is a unit matrix andλ a (a = 1, 2, 3) are the Pauli matrices. First, we consider the mass spectrum of dark pions by use of the LSM. In the absence of the source term H, dark pions are massless Nambu-Goldstone bosons. H tilts the LSM potential and gives the masses of dark pions. In this article, we consider the case j 0 dominates the source term with a small perturbation j 3 that describes isospin violation in dark quarks, while we also discuss other cases in appendix A. In this case, σ 0 obtains its VEV of order of the dynamical scale √ N C Λ χSB /4π, while σ 3 has its VEV due to a tilt by a source term j 3 σ 3 , and thus the VEV σ 3 is proportional to the source j 3 Λ −2 χSB . The spectrum of pseudoscalars is as follows.
Here, f π ≡ σ 0 (4π) −1 √ N C Λ χSB corresponds to the pion decay constant, and δ ≡ σ 3 . π 0 mass is dominated by the anomaly term, which corresponds to the η meson in the SM, while the mass of other pseudoscalars is proportional to j 0 at the leading order. The mass difference within the isospin multiplet is proportional to square of the dark quark mass difference, (M 1 − M 2 ) 2 (see e.g., [101,102]).
Besides the quark mass difference, the U(1) D interaction explicitly violates the isospin symmetry and provides the pion mass difference. The effective Lagrangian up to O(e 2 ) is given by where c π q is an O(1) dimensionless parameter, and α = e 2 /4π. This Lagrangian is invariant under the global flavor symmetry as the charge matrices q L and q R are spurions. We -8 -

JHEP03(2022)176
determine the coefficients to give a contribution only to the charged pion mass when the charge matrices are in charge basis (q L = q R =q).
Here, the charges of dark quarks determine the charge matrix in the charge basis,q. 2 The mass differences among dark pions originate from U(1) D interaction and dark quark mass difference, and each of them is numerically given by (2.21) , and large-N C scaling of other parameters. m π denotes the mass of π 3 . As with the pion mass difference in the SM (see ref. [103]), the U(1) D correction dominates the dark pion mass difference unless dark quark mass difference is quite large.
After diagonalizing the dark pion mass terms including the U(1) D correction, the dark photon couples to dark pions.
(2.22) Here, we use primed notation for the hadronic field in mass basis after taking into account the U(1) D correction. An additional angle β V arises from diagonalizing the mass terms including the U(1) D correction in addition to θ V . We discuss the dark pion mass spectrum in detail in section B.
Next, we consider the spectrum of light dark nucleons. The Lagrangian density of the dark nucleons is A typical scale of the dark nucleon mass is determined by the dynamical scale, gf π /2, and the NDA with the large-N C scaling of the coupling g is g 4π √ N C . We take into account two mass terms when we consider the mass difference of the dark nucleons: one comes from the VEV of the LSM field, and another comes from the dark-quark mass matrix itself which is proportional to an O(1) dimensionless parameter ξ. The mass difference that 2 The charge matrix of dark quarks is determined by T 3 + B/2 with isospin T 3 and baryon number B of dark quarks. The baryon number is zero for Φ, and the charge matrix of Φ can be determined only by the isospin. However, since Φ is identified with Q Q , the charge matrix of Φ can also be determined byq, the charge matrix of Q . Since we can add any matrix proportional to unit matrix to the charge matrix in charge basis, both definitions of the charge matrix do not change the results (except for anomaly induced interactions). It is easy to compute the spectrum by use of the charge defined by T 3 + 1/2 = diag(1 , 0).
comes from the VEV of the LSM field is proportional to j 3 , which originates from the mass difference of the dark quarks. When j 0 is dominated, the dark nucleon masses are given by The lightest dark nucleon is the candidate of the composite ADM. As shown in refs. [67,85,104], if the asymmetry is fully shared between the dark and SM sectors, the DM mass is to be 8.5/n g GeV with n g being the number of generations of dark quarks. When we consider mirror-sector scenario to explain the coincidence of the dynamical scales in the visible and the dark sectors [105], n g should be integer. As long as we consider vector-like QCD in the dark sector, n g can be half-integer. The one-loop beta function of the dark QCD coupling is proportional to 11−4n g /3. The dark QCD is asymptotically free as long as n g ≤ 8, and this gives a rough lower bound on the ADM mass to be larger than 1.1 GeV. On the other hand, the large number of generations causes the Landau pole of the dark QED coupling. The one-loop beta function of the dark QED coupling is proportional to 20n g /9 + 1/3. The dark QED coupling does not diverge up to the Planck scale when we take the dark QED coupling to be α 0.06/(n g + 0.15) at 1 GeV for n g generations. In the following, we assume N 2 is the lightest nucleon, and we set m N 2 = 8.5/n g GeV.
In addition to the mass difference from δ and dark quark mass difference, we can add the U(1) D correction as with the dark pions. The U(1) D correction is typically α f π /4π in the limit of small θ V and δ . (2.25) As for the dark nucleons, we ignore the U(1) D correction in the following since ∆ QED N is smaller than ∆ iso N . The dark hadrons with different U(1) D charges can mix with each other since U(1) D is broken in our model. In this case, the heavier dark hadrons mainly decay into the lighter ones and (off-shell) dark photon. The charge matrices in the mass basis are important for determining the decay rate. As we assume the quark mass matrix M Q to be hermitian in this article, the diagonalizing matrices are identified, U L = U R = U V . The dark hadrons in the charge basis are obtained by the same unitary matrix U V : This rotation also gives the relation between the charge matrices in the mass basis and in the charge basis as follows: (2.28) We ignore the phase α V in the following. Since we assume that U(1) D is spontaneously broken by the dark Higgs φ D , the dark photon mass is not related to the mass spectrum of dark hadrons. We take the mass of dark photons as a free parameter in this study. The dark photon mass is related with the dark pion mass, for instance in a chiral setup of the dark confining sector [58,106,107].

Dark hadron mass spectrum and lifetime
We have discussed the generic mass spectrum and interactions in the composite ADM model with a dark photon in the previous section. The search strategies for the long-lived particles that are produced at the collider experiments depend on their lifetime.
When the produced particles are quite long-lived or are stable in the dark sector, they leave missing signals or signals from scattering with the detector material. In the composite ADM models, the dark proton which has a U(1) D charge scatters with the SM particle via dark photon, and then the scattering gives the recoil energy of nuclei in the detector material [108]. When the dark neutron that does not have a U(1) D charge is stable, it escapes from the detectors, and then this process has been and will continue to be explored by searches for invisible decay at NA64 [109,110], LDMX [22,23], and Belle-II [24,25].
The dark pions are lightest among the dark hadrons. When the dark pions are lighter than dark photons, the dark pions decay into visible particles with off-shell dark photons. In this case, the dark pions are quite long-lived, and hence their decay leads to problematic energy injection in electromagnetic channels at the late time, in particular, when their lifetime is longer than 10 4 s [111]. As for the dark photon decay, we have different phenomenology when a decay mode A → π π opens. The decay has been investigated by the invisible dark photon decay searches.
On the other hand, in this study, we focus only on the visible decay signals from dark hadrons. The dark hadrons should have decay length of O(1) m in order to leave visible signals at the detector. To begin with, we discuss the lifetime of the dark hadrons, and determine the mass spectrum that is the most relevant to the visible decay searches. Then, we fix the nucleon mass spectrum from the direct detection experiment for the dark matter. Dark photons are lightest in the dark sector, and thus we discuss the dark photon searches in the end of this section. We will discuss the production of dark hadrons at the collider experiments and the fixed-target experiments in section 4.

Decay and transition of dark hadrons
When the mass difference of the dark nucleons is smaller than the dark photon mass m A , the heavier dark nucleon N 1 mainly decays to the lighter nucleon N 2 with the SM fermions.

JHEP03(2022)176
The differential decay rate is given by Here, s ff ≡ (p f + pf ) 2 . Γ N 1 →A N 2 is the decay width into on-shell dark photon with the dark photon mass replaced with m A → √ s ff . Γ A denotes the decay width of the dark photon with the mass of √ s ff , which is given by eq. (2.4). We include the hadronic final state when the fictitious mass √ s ff exceeds the threshold. The momentum transfer is much smaller than the dark photon mass, s ff m 2 A , when the decay mode into on-shell dark photon does not open. In this case, in the small limit of s ff and the mass difference, denoted by where α = e 2 /4π is the fine structure constants of U(1) D . Once the mass of the SM fermions is neglected, the decay rate is approximately given by Here, α = e 2 /4π is the fine structure constants of the electromagnetism and Q f is the electromagnetic charge of f . The approximate decay length is given by (3.4) The lifetime of the heavier state is much smaller than 1 s, and then the heavier states that are produced in the early universe do not lead to the problematic energy injection to the electromagnetic channels at the late time. Therefore, their decay is free from the cosmological constraints, and we focus on the visible decay searches and will discuss future sensitivities to the decay in section 4.1.
Similarly, the heavier dark pion (with mass m π 1 ) can decay into the lighter one (with mass m π ) by emitting the off-shell dark photon. Once we neglect mass of final state fermions, the decay rate through the dark pion transition is given by Here, ∆ π = (m π 1 −m π )/m π and 2θ V +β V denotes the mixing angle between charged and neutral dark pions, which we discuss in section B. However, since dark pions are lighter than dark nucleons, decay length of pions through transition tends to be much longer than that of nucleons. Unless pion mass is close to nucleon mass, the dark pion decay via the pion transition is expected to be beyond the sensitivity of the visible decay searches. Therefore, we do not include the dark pion transition in this study. A dark-neutral pion can couple to dark photons through the anomaly-induced interaction: (3.6) Here, f π is a decay constant of the dark pion. The dark pion mainly decays into (off-shell) dark photons via this interaction: as illustrated in figure 1, decay into two on-shell dark photons as 2m A < m π , decay into a dark photon and a pair of the SM fermions as m A < m π < 2m A , and decay only into the SM fermions as m π < m A . The decay rates to two dark photons and to a dark photon and a pair of the SM fermions are, respectively, given by Here, the decay rate for π 3 → A ff is computed neglecting the fermion mass, and F (x) is a function originated from the three-body phase space integral: (3.9) This function vanishes as x → 1. The dark pions promptly decay into to dark photons when 2m A < m π . On the other hand, the dark pions are long-lived when m A < m π < 2m A due to the three-body phase space and the kinetic mixing . When we set m π = 1.8m A , a proper decay length of the dark-neutral pions is given by (3.10) When the dark photon is heavier than the dark pion, all final states of the decay process are the SM fermions. As shown in figure 1, we have two processes only with the SM fermion final states: one is the four-body tree-level decay, and another is the loopinduced two-body decay. We use the four-body decay rate in the massless fermion limit Figure 2. Decay lengths of dark-neutral pions: the kinetic mixing is assumed to be = 10 −2 (solid), = 10 −3 (dashed), and = 10 −4 (dotted). We take α = 0.05 , m π = 0.5 GeV , f π = 0.8 GeV. We include all possible final state for π 3 → A ff , but include only leptonic final states for π 3 → ff and π 3 → ff f f . and the two-body decay rate that are computed in ref. [112] as dark η meson decay. We take into account the difference of the anomaly coefficients as the dark pion decay.
Due to the extra suppression from kinetic mixing 4 and the four-body phase space or the loop suppression factor, the decay rate without an on-shell dark photon final state is very tiny. Hence, the visible decay searches are not promising when the dark pion is lightest in the dark sector. We show the decay length of the U(1) D -neutral dark pion in figure 2 with the dark pion mass of 0.5 GeV and the decay constant f π = 0.8 GeV. The different line-types in the figure correspond to different values of kinetic mixing : = 10 −2 (solid), = 10 −3 (dashed), and = 10 −4 (dotted). Each lines is composed of three parts: prompt decay via π 3 → A A for 2m A < m π , long-lived decay via π 3 → A ff for m A < m π < 2m A , and very long-lived decay via four-body/loop-induced decay. As shown in figure 2, dark pions would have the typical decay length of O(10 −2 )-O(10 2 ) m via three-body decay when m π is close to 2m A . Therefore, the visible decay searches have a potential to explore the dark pion decay in this mass range.
The darkly charged pions can also decay through the pion mixing and chiral anomaly. When the darkly charged pions decay into the SM particles, their lifetime is longer than -14 -

JHEP03(2022)176
dark-neutral pions due to the extra factor of the pion mixing. When the kinetic mixing is of order of O(10 −3 ), π ± → A A will give the visible signal since the dark photon promptly decays into the visible particles in this range of the kinetic mixing. To get the decay length of O(1) m, the pion mixing angle should be or O(10 −4 ) or below. Meanwhile, we assume the dark nucleon mixing θ V to be O(10 −2 ) in order that the decay length of the dark nucleon transition is of O(1) m. The dark pion mixing is expected to be the same order of magnitude as the dark nucleon mixing, and hence we focus only on the decay of dark-neutral pions in this study.
In this article, we focus on the visible decay searches, namely we focus on m A ≤ m π ≤ 2m A . We will discuss the sensitivity of the searches for the dark pion decay in section 4.2.

Direct detection
The DM has a coupling to the dark photon through the nucleon mixing as shown in eq. (2.28), and hence the DM direct detection experiments put constraint on the kinetic mixing . The U(1) D charge of the lightest dark nucleon is sin 2 θ V , where θ V is the nucleon mixing angle. The predicted ADM mass is m N 2 = 8.5/n g GeV, and in particular we consider m N 2 = 8.5 , 4.3 , 2.1, and 1.1 GeV (corresponding to n g = 1 , 2 , 4 , and 8 , respectively) in the following. The scattering of the ADM off nuclei leaves recoil energy that is close to the energy threshold of the direct detection experiments with liquid noble gas. For this reason, we use the constraints given by the direct detection experiments with unconventional methods, in particular for the light ADM.
The current bounds on the DM-nucleon cross section is the following. The direct detection experiments using liquid xenon constraint on the DM-nucleon cross section get weaker for the DM lighter than 10 GeV. PandaX-II experiment puts constraints on the DMnucleon cross section both via a heavy mediator (contact interaction of DM and nucleons) and via a light mediator [113,114]. When the DM mass is 8.5 GeV and the dark photon mass is 1.0 GeV, the upper bound on the cross section is given by σ ≤ 6 × 10 −45 cm 2 . Comparable constraints are obtained from XENON1T (S2-only technique, using ionization signals) [115] and LUX [116], and one order of magnitude severer constraints are obtained from XENON1T (with the conventional method) [117]. We use the constraints on the spinindependent cross section for the DM mass of 4.3 GeV and 2.1 GeV by DarkSide-50 [118], which uses the liquid argon. The corresponding bounds on the DM-nucleon cross section are given by σ ≤ 2 × 10 −42 cm 2 for the DM mass 4.3 GeV and σ ≤ 10 −41 cm 2 for the DM mass 2.1 GeV. Even though we use the constraints by DarkSide-50, as for the DM mass of 4.3 GeV, XENON1T has put a constraint one order of magnitude severer than that of DarkSide-50 [115]. Direct detection experiments with solid-state detectors have placed constraints on the DM-nucleon cross section even for the DM mass of GeV. The CRESST-III experiment has already set an upper bound on the cross section to be less than 5 × 10 −39 cm 2 for 1 GeV DM mass [119], which is about six orders of magnitude larger than the constraint for the ADM mass of 8.5 GeV. A light dark matter search at XENON1T (with the use of Migdal effect [120][121][122]) has reported, and the constraint is comparable to the CRESST-III experiment [123].

JHEP03(2022)176
When the dark protons are the lightest dark nucleon, the bound on the kinetic mixing is given by ≤ 1.4 × 10 −7 (m A /1 GeV) 2 as the DM mass m N 2 = 8.5 GeV and the fine structure constant α = α [67]. As we will see in the next subsection, this bound is far below the future sensitivity to the visible decay searches of dark photons. In contrast, direct detection bounds get much weaker when dark-neutral baryons, i.e., dark neutrons, are the lightest nucleons since dark neutrons do not have an electric coupling to dark photons. 3 Indeed, a dark-neutral baryon is essential in composite ADM models, in order to share the generated asymmetry between two sectors via the portal interactions eq. (2.2).
In our study, the lightest dark nucleon is mainly composed of the dark neutron, and mixes with the dark proton. The bounds from the direct detection experiments are replaced by the bound on the kinetic mixing and the nucleon mixing as follows. (3.12) Here, σ bound denotes the upper bound on the nucleon-DM scattering cross section given by PandaX-II, DarkSide-50, and CRESST-III as shown above. The upper bound on the kinetic mixing depends also on the detector material as A/Z, which is taken to be that in PandaX-II, namely xenon. The upper bound on the kinetic mixing is relaxed thanks to the small mixing θ V . As for the ADM with the mass of 8.5 GeV and small m A , we take into account the constraints on the scattering cross section with a light mediator given by PandaX-II [113,114].

Dark photon searches
The search strategy of the dark photon visible decay changes depending on lifetime. Various prompt decay searches have constrained the dark photon parameters. For the dark photon mass above GeV, dark photons can be produced via Drell-Yan production and heavy rare meson decay in collider experiments, and hence collider experiments put constraints on the kinetic mixing from the absence of the signals. LHCb has performed the inclusive searches for the dark photon (e.g., pp → XA → Xµ + µ − with any final states X), and put the upper bound on the kinetic mixing of 10 −3 from prompt-like decay search [27,125,126]. At electron colliders, dark photons are produced by meson decay (e.g., φ → ηA ) and radiative return process (e + e − → γA ). The absence of a deviation from the SM background has been already reported at BaBar [8] and KLOE [9,[80][81][82] and has set upper bounds on the kinetic mixing in a broad dark photon mass range from MeV to GeV.
There have been several bounds on dark photon parameters from astrophysical observations. New light particles with weak coupling to the SM particles can be produced in the hot supernova core via mainly bremsstrahlung, and can be radiated. This process provides a new cooling process of supernova, and leads to inconsistency with the observation of neutrinos from supernova, SN1987A [127][128][129][130][131]. A typical decay length of the dark photon that is constrained by SN1987A is about 10 km corresponding to the size of the supernova core. The dominant bound arises from the supernova cooling for dark photons with the mass below 100 MeV and with the kinetic mixing of O(10 −9 -10 −7 ), and therefore the parameter space constrained by SN1987A is beyond our interests.
Fixed-target experiments have the sensitivity to the parameter space between constraints by the prompt decay search in collider experiments and the SN1987A observation. Electron beam-dump experiments, which have the beam dump area without any tracking system of daughter particles, put constraints on the parameter space above 10 −7 and below a few hundred MeV. The dark photons are produced via Bremsstrahlung at the electron beam-dump experiments. The strong constraints have mainly come from the E137 [12] and E141 [132] experiments at SLAC, the Fermilab E774 experiment [133], Orsay [134], and KEK [135]. Though they have mainly searched for axion-like particles and light Higgs bosons, the constraints on the dark photon models are obtained by reanalyzing these constraints [4,7]. The electron fixed-target (non beam-dump) experiments have a sensitivity between those of prompt decay experiments and beam-dump experiments since their detectors are not located very far from targets. There have been the electron-beam fixed target experiments searching for visible final states: A1/MAMI [136], APEX [137,138], and HPS [16,139,140]. Further running of APEX and HPS is planned in the near future.
Similarly to the electron beam-dump experiments, data from previous proton beamdump experiments have been also reanalyzed to constrain the dark photon visible decay: CHARM [141,142], LSND [5,143], and U70/νCal [144,145]. The use of proton beam provides the enhancement of dark photon production, due to high penetration rate of protons on target materials and enhanced production mechanism from meson decays besides Bremsstrahlung processes. The on-going and future fixed-target experiments with proton beam, SeaQuest [20,21] and SHiP [17,18], have been proposed. Since the decay volumes of the proposed experiments are shorter than those of the past proton beam-dump experiments, the proposed experiments have the sensitivity to dark photons with short lifetime.

Dark hadron searches
In this section, we discuss visible decay searches of dark hadrons at the lifetime frontier. In particular, we consider the dark nucleon transition and the three-body decay process of the dark-neutral pion. We assume that the dark hadrons made of one of n g generations dark quarks are predominantly produced at the LHC and the fixed-target experiments. The dark hadrons are produced via off-shell dark photons at the LHC and the fixed-target experiments since the dark photons are lightest in the dark sector. The differential cross section for the production of dark hadrons is written as [108] dσ

JHEP03(2022)176
Here, x collectively denotes the variables for the dark photon production, and m * A is the fictitious dark photon mass corresponding to the invariant mass of final states. σ A denotes the cross section of on-shell dark photon with the dark photon mass m * A being the fictitious mass, and Γ A is the partial decay width of A with the mass of m * A . From this formula, we evaluate the produced number of dark hadrons via the off-shell dark photons when the invariant mass is above the dark photon mass as follows. (4.2) Here, N A denotes the number of produced dark photons. The visible signal would be detected when the long-lived particles decays into the SM particles inside the detector volume or generically in a specific region. The number of signals is roughly evaluated as follows: Here, we assume the branching fraction of the long-lived particles to the visible particles to be unity. (eff) denotes the total efficiency when the geometric detector acceptance has weak dependence on the location of the decay: A eff is the acceptance of detector, r min and r max are, respectively, distances of the detector starting and end points from the interaction point, and d ≡ cτ pm −1 is the boosted decay length of the long-lived particle with the boost factor pm −1 , and . . . denotes the ensemble average. Though the efficiency depends on m * A , the fictitious mass of off-shell dark photon, we factorize the efficiency (eff) at m * A where the signal number N signal is dominated. The acceptance A eff includes the geometric acceptance A geo and the depletion of the produced numbers at a certain momentum differs from typical momentum A prod : A eff = A geo A prod . Since the detectors are located away from the interaction points at collider or the beam-dump areas, the angle of the emitted particle with respect to the beam-axis is essential to detect the visible signal. The angle of the emitted particle is determined by the decay of off-shell dark photon into the longlived particles since the emission of the off-shell dark photon is collinear to the beam axis in the following. Meanwhile, the energy threshold of detectors is imposed for some reasons discussed later, for instance to remove low-energy backgrounds. When the longlived particles predominantly coming to a detector have the energy below the threshold, most of long-lived particles is vetoed, and the number of long-lived particles leaving visible signal decreases. We will discuss the acceptance that we use in this study later in detail.
In the following, we consider the LHC experiment and the fixed-target experiments. The dark photons are mainly produced through the Drell-Yan process in the former and the Bremsstrahlung process in the latter. We use the produced number of dark photons through the Drell-Yan process at the LHC [146] and through the Bremsstrahlung process at the electron beam fixed-target experiment, E137 [108], and at the proton beam fixedtarget experiment, SeaQuest [21]. The produced number N A via the Drell-Yan process and the electron Bremsstrahlung process is proportional to 2 m −2 A , while N A via the proton -18 -

JHEP03(2022)176
Bremsstrahlung does not depend on the dark photon mass for m A 500 MeV. As for the proton Bremsstrahlung, N A is enhanced due to the mixing of A and the SM ρ-meson near m A m ρ [147], and N A is sharply suppressed due to the form-factor for m A above 1 GeV, where the internal structure of nucleons is essential. 4 When the long-lived particle is elementary, the decay rate Γ A is proportional to m A , and therefore the integral given by eq. (4.2) is dominated by the infrared (IR) contribution for the electron Bremsstrahlung. As discussed in ref. [108] for the electron Bremsstrahlung, through the off-shell dark photon, the dark-sector particles are dominantly produced at m * A near the kinematical threshold 2m χ with m χ being the dark-sector particle mass: . On the other hand, when we consider composite particles, the decay rate can have a higher power of m A , and then the UV contribution (e.g., dark dynamical scale) can dominate the produced number N . We will discuss the UV dominance in the dark hadron production later. As for the proton Bremsstrahlung, N is predominantly determined by the SM ρ-meson resonance even when we consider composite particles due to the enhancement of N A at the ρ-meson resonance.

Dark nucleon searches at the LHC
We discuss the visible searches at the LHC lifetime frontier for the transition of the dark nucleon. The off-shell dark photon is produced via the Drell-Yan process, and then N 1 , the heavier dark nucleon, is produced via the dark hadronization when a sufficient energy through the off-shell dark photon is injected to the dark sector. Various experiments have been proposed to explore long-lived particles with the decay length of O(1)-O(10 2 ) m at the LHC: MATHUSLA [29], FASER [30], and CODEX-b [31]. First, we summarize the LHC lifetime frontier.
The MATHUSLA experiment has been proposed [29] to search for long-lived particles at the LHC. Its detector would be placed on the surface ∼ 100 m above the LHC beam line and ∼ 100 m downstream from the interaction point, and would take data during the HL-LHC. The MATHUSLA detector consists of on-top five tracking layers with timing resolution of 1 ns, air-filled decay volume, and floor-detector reducing background muons from the CMS interactions point. Although, different decay volumes, so called MATHUSLA50, MATHUSLA100, and MATHUSLA200, have been considered in the MATHUSLA proposal [29,149,150], we focus only on MATHUSLA200 with decay volume of 200 m×200 m×20 m in this article. At the MATHUSLA experiment, various backgrounds are expected to be reduced by the following. The rock is present between the LHC ring and the detector location, good time resolution of the tracking layer differentiates upward going (long-lived particle) signals from downward going (cosmic ray) signals. Therefore, the searches for the long-lived particles are assumed to be background free at the MATHUSLA.
The FASER experiment has been proposed [30] to search for long-lived particles in the very forward direction of the ATLAS interaction point. Its detector is located in the very forward region along the beam line of the LHC, and is positioned ∼ 480 m away JHEP03(2022)176 from the interaction point, where the LHC tunnel just starts to curve. The detector has a cylindrical shape, and its longitudinal axis coincides with the LHC beam collision axis. Various versions of the FASER experiment have been proposed: FASER [151,152] and FASER2 [153]. The former collects data during LHC Run 3 and its detector have the cylindrical shape with 1.5 m depth and 10 cm radius, while the latter takes data during HL-LHC with the cylindrical shape of 5 m depth and 1 m radius. As with the MATHUSLA experiment, background events are expected to be reduced by the following. There are the 100 m of rock and concrete shield of the detector from the ATLAS interaction point, and a charged particle veto in front of the detector. We focus on the future sensitivity of the dark nucleon at the FASER2 experiment.
Similarly to the MATHUSLA experiment, it has been proposed to put a new detector near the LHCb interaction point, which is named as the CODEX-b experiment [31]. The detector would be placed about 25 m away from the interaction point, with a 10×10×10 m 3 volume. The background events are expected to be reduced by the existing concrete wall of ∼ 3 m, an additional lead shield of ∼ 4.5 m, and a charged-particle veto inside the lead shield [154]. In comparison with the MATHUSLA detector, CODEX-b has a smaller detector volume, the detector location is closer to the interaction point, and the integrated luminosity at LHCb is smaller than that at ATLAS. Due to the geometry and the expected luminosity of CODEX-b, the future sensitivity of CODEX-b is comparable with but lower than that of MATHUSLA in the inelastic dark matter model [146]. Thus, we do not include the future sensitivity of CODEX-b in the following analysis.
The searches for long-lived particles at the LHC have been studied in the context of the inelastic dark matter model in refs. [26,146]. We map their results into the dark nucleon transition (the dark-neutral pion decay in the next subsection) by taking into account important differences. To this end, we outline their setup to compare with our model. In these studies, U(1) D gauge invariance allows only a Dirac mass term of two Weyl fermions, and Majorana masses arise from the spontaneous U(1) D breaking. In the mass basis, the lightest state χ 1 is the DM and inelastically couples to the heavier state χ 2 with dark photons. Their mass splitting is controlled by the nearly degenerated Majorana masses, and the mass splitting is small when the Majorana masses are smaller than the Dirac mass. Dark photons are assumed to be heavier than DM particles; in particular they assume m A = 3m χ where m χ is the DM mass. χ 2 is produced via dark photon decay, and hence the produced number of χ 2 denoted by N 0 χ is equal to that of dark photons: N 0 χ = N A . Then, the produced χ 2 decays into χ 1 and the SM fermions via off-shell dark photons.
In our study, we consider that dark photons are lighter than dark nucleons as discussed before. In this case, we have an additional factor for production rate of N 1 via the dark hadronization as follows: 5

JHEP03(2022)176
Here, N N 1 denotes the number of N 1 and N 1 produced in collider experiments. The factor 5/3 originates from the sum of the charge squared of dark quarks. We take m * A = √ 10m N 1 by analogy with the off-shell production of the elementary dark particles [108]. This may be a crude approximation since it is not clear if the contribution near the kinematic threshold dominates the production even for the composite dark particle. n N 1 denotes the multiplicity of N 1 for the production through the dark photon mediation. For simplicity, we take a similar value, n N 1 = 0.04, to the multiplicity of protons by the e + e − annihilation for the center of mass energy just above J/ψ threshold [155]. n N 1 depends on the choice of m * A since the multiplicity depends on the injected energy. We discuss how the change of the multiplicity affects the sensitivity of the visible decay searches at the LHC lifetime frontier in section D.
The sensitivity range of the kinetic mixing has a maximum value, max , and a minimum value, min , with other parameters being fixed. Due to detector designs of MATHUSLA and FASER, the visible decay products are vetoed as background events if N 1 decays before the decay volume area. This means that there is a lower limit on the decay length of N 1 , and below this limit N 1 cannot leave any visible signal at the detector even when N 1 is maximally boosted. It leads to the maximum value max . The maximum value max is evaluated by use of eq. (4.3) for the boosted decay length d r min < r max : Here, p * max (subscript corresponds to max , and does not necessarily coincide with the maximum energy) is given by the maximum momentum, which we will discuss later. Instead of directly computing N N 1 and A eff , we consider a formula similar to eq. (4.6) in the inelastic DM, with max and other parameters given in ref. [146]. By using the similar formula in the inelastic DM, we can rewrite eq. (4.6) as follows: Here, d 0 (p max ) is the boosted decay length of χ 2 computed with the reference parameter set in the inelastic DM that corresponds to a reference point on the upper boundary max of sensitivity area in a m 1 -plot read out from ref. [146] and with the momentum p max same as p * max . In this formula, N N 1 appears as a ratio to N 0 χ , and therefore we only have to care about the parameter dependence of N A in eq. (4.5). In particular, N A is proportional to 2 m −4 A for FASER (A only in the forward region) and 2 m −2 A for MATHUSLA. A eff and A 0 eff denote the efficiencies in the composite ADM model and in the inelastic DM model, respectively. The efficiency also appears as a ratio to A 0 eff , and hence we only have to care about the model parameter dependence of the efficiency. We will discuss the efficiency in the end of this subsection in detail. We demonstrate the validity of our approximation in appendix C by applying the formula to the other parameter sets in the inelastic DM model.
Concerning p * max , we take p * max = 7 TeV for FASER and p * max = 100 GeV for MATH-USLA. We expect that p * max hardly depends on the model parameters. Dark nucleons are assumed to be mainly produced by the Drell-Yan process at the LHC lifetime frontier. Due to the parton distribution function, the production cross section through the -21 -

JHEP03(2022)176
Drell-Yan process sharply decreases above about 100 GeV (e.g., see ref. [156]). This implies that it would be hard to produce the long-lived particles with the transverse momentum above 100 GeV at LHC collision. Meanwhile, in the forward direction to the beam line, the produced particles can carry the momentum of the beam energy. Therefore, we take p * max = 100 GeV (7 TeV) for MATHUSLA (FASER). The long-lived particles hardly decay inside the detector volume when the boosted decay length is larger than r max : r min < r max d. This determines another boundary of the sensitivity plots. In this case, min satisfies where p * min (subscript corresponds to min , and does not necessarily coincide with the minimum energy) denotes momentum of N 1 leaving signals, which we will discuss later. Similarly to max , we estimate min by considering a formula similar to eq. (4.8) in the inelastic DM model instead of directly calculating N N 1 and A eff . By using the similar formula in the inelastic DM, we can rewrite eq. (4.8) as follows: Here, d 0 (p min ) is the boosted decay length of χ 2 computed with the parameters in the inelastic DM that corresponds to from ref. [146] and with the same expression of momentum p min as p * min (with proper replacement of model parameters). Only the ratio N N 1 /N 0 χ appears again, and hence we do not care about the overall factor to the produced number, and we estimate it by N A and its parameter dependence. A eff appears only as a ratio to A 0 eff . The validity of our approximation is demonstrated in appendix C by applying our formula to the inelastic DM model. This approximation formula determines min , and hence it seems to prefer to take p * min to be the minimum value with which the particles leave the visible signal above the energy threshold at the detector, p thr , since it leads to the least boosted particles. However, the signal number of N 1 with a momentum is suppressed when the momentum of N 1 deviates from a typical momentum of N 1 leaving signals. N 1 leaving signals typically has the momentum p * geo (m * 2 A − 4m 2 N 1 ) 1/2 /2θ det where θ det is the angle of N 1 with respect to the beam axis. We basically take p * min = p * geo to make the acceptance maximum as much as possible, but p * min does not coincide with p * geo when the particle with the typical momentum cannot leave a visible signal due to the energy threshold (i.e., p * geo < p thr ). We will discuss the energy threshold for the visible decay products at each experiment later. Therefore, we take p * min = min(p * max , max(p thr , p * geo )) unless p thr exceeds p * max where nothing is detected. The produced dark nucleons have to be energetic since decay products should be sufficiently energetic because of energy threshold at the detectors. Heavier dark nucleons decay into a lighter dark nucleon and visible particles via small mixing angle. The final state dark nucleon takes over most of the energy that initial dark nucleon has, and then the SM particles can be less energetic. Hence, the energy thresholds play a significant role to determine the sensitivity. In this study, we take the energy threshold at MATHUSLA and FASER as follows. (4.10) In ref. [157], the final state energy thresholds at MATHUSLA are discussed. The thresholds range between E thr min 200 MeV-1 GeV in order to detect displaced vertices efficiently at MATHUSLA, and then we take E thr min 600 MeV per each track to follow the analysis by ref. [146]. Meanwhile, the background events at FASER have been well discussed in ref. [151]. Following the literature, we assume E thr min 100 GeV for the total energy deposition in order to reduce the trigger rate and to remove low-energy backgrounds at FASER [151,153]. The corresponding momentum of N 1 , which is denoted by p thr , is given by (4.11) Here, we multiply the factor 2 for MATHUSLA to make all the SM decay products sufficiently energetic since the energy threshold is imposed on each charged track. Since the threshold value p thr exceeds the maximal value, the FASER experiment decreases its sensitivity for the dark nucleon search if ∆ N 0.03.
In the case of MATHUSLA, for the fictitious dark photon m * A 10 GeV, the typical momentum that is close to the mass of dark photon is smaller than the threshold. In this case, the signal number of N 1 decrease and depends p thr /p * geo after the final state phase space integration. We empirically find the efficiency A eff related with the decrease of the signal number from [146] as follows.
A prod (p * ) is normalized to be unity as the long-live lived particle has the typical momentum. Though we take the MATHUSLA case as an example, we take into account the depletion of signals as A prod (p * ) even for the FASER case when the momentum deviates from p * geo . As far as we take p * min = p * geo , the produced particle comes in the detector. However, it may not come in the detector when p * min deviates from p * geo . We incorporate the angle acceptance as follows.
and θ det is determined by each experiment: θ det = 0.5 for the MATHUSLA and θ det = 2 × 10 −3 for FASER. The total efficiency for the LHC lifetime frontier is defined by We summarize the parameters featuring the experimental setups of MATHUSLA and FASER in table 2. As for p * min in this table, we show just typical values: p thr fro MATH-USLA and p * geo for FASER. In this table, we also show the parameters for the dark pion searches, which will be discussed in the next subsection.

Dark pion searches
We consider the dark-neutral pion decay with a specific mass spectrum m A < m π < 2m A in this subsection. The dark-neutral pions can be produced via various mechanisms: 1) transition from the darkly-charged pions, 2) dark hadronization, and 3) via off-shell dark photon with the Wess-Zumino-Witten (WZW) term [158,159]. The dark-neutral pions can be produced via the transition from darkly-charged pions that are produced from offshell dark photon as in the dark nucleon production. However, the transition rate from the darkly-charged pions is considerably tiny, and then the produced darkly-charged pions escape from the detectors without any signals. In the SM, a hadronic channel π + π − π 0 is opened with injected energy above the ω-meson threshold [89,160,161]. In a similar manner, dark-neutral pions are produced via the hadronization when the injected energy is larger than the dark dynamical scale. We use the production via the dark hadronization for the LHC lifetime frontier, and use it even for fixed-target experiments when the injected energy is larger than the dark dynamical scale. When the injected energy is less than the dark dynamical scale, we utilize the production via off-shell dark photon. We discuss each dark-pion production process in detail. As for the dark hadronization, we use an approximation formula for the dark pion production, N π n π α π (4.14) Here, n π denotes the multiplicity of dark-neutral pions for the production through the dark photon mediation. The factor 5/3 originates from the sum of the charge squared of dark quarks. We take a similar value, n π = 2.0, to the multiplicity of pions by the e + e − annihilation for the center of mass energy just above J/ψ threshold [155]. Once the injected energy to the dark photon exceeds the dynamical scale, the dark hadrons are expected to be produced through the hadronization. In particular, the dark-neutral pions are expected to be produced above dark-isospin-singlet vector resonance (corresponding to ω-meson in JHEP03(2022)176 the SM), and therefore we take the fictitious dark photon mass to be m * A = m N 2 . We discuss how the change of the multiplicity affects the future sensitivity to the dark pion decay at the LHC lifetime frontier in section D.
We introduce the relevant fixed-target experiments before discussing the production of dark pions at each visible decay search. As discussed in section 3.3, there exist data from previous fixed-target experiments that have excluded visible decay of sub-GeV dark photons with 10 −7 -10 −5 . The same data can be applicable for constraints on the three-body decay of the dark-sector particles when the final states include visible particles. The decay length of three-body decay is enhanced compared to that of two-body decay due to the final-state phase space, and hence a search for three-body decay favors the long-baseline fixed-target experiments. Among the fixed-target experiments with electron beam, the E137 experiment is sensitive to longer lifetime than others since it has a very long natural shielding (hill of 179 m) and huge decay volume before detector (open air region of 204 m). The E137 experiment utilizes a 20 GeV electron beam, and the beam is dumped in a wateraluminum target. No signal events with the energy deposit above 3 GeV were observed at E137, and therefore the expected signals from the dark-sector particles should be less than a few. CHARM and U70/νCal that are existing fixed-target experiments with proton beams have constrained the similar parameter space as E137. The visible search at E137 has put a constraint on dark photon parameters [4], and hence we map the results of the existing constraint for the visible decay of sub-GeV dark photon into the visible decay of the dark pions. In the literature, the on-shell dark photon is produced via the electron Bremsstrahlung, and the produced dark photon decays into the SM particles via kinetic mixing.
The SeaQuest experiment is currently running at Fermilab with the 120 GeV proton beam [20], and is originally designed to measure antiquark structure of nucleons via Drell-Yan dimuon production. A magnetized iron block of 5 m is placed just downstream from a target to sweep away the soft SM radiation, and the detector is composed of a 3 m magnet and four tracking stations, which allow us to reconstruct decay vertex and momentum precisely. The first tracking station is located 1 m downstream from the beam-dump, and hence the decay volume of the SeaQuest is in the range of r min = 5 m and r max = 6 m. There is a plan to upgrade SeaQuest by installing an electromagnetic calorimeter that was utilized in PHENIX experiment at Brookhaven National Laboratory, which is named as "DarkQuest". This upgrade allows to detect electrons as decay products of dark-sector particles. Since the decay volume of SeaQuest is shorter than the other fixedtarget experiments, the SeaQuest is sensitive to the dark photon parameter for shorter decay length. The visible signal search at SeaQuest has been studied in the context of the strongly interacting massive particles (SIMP) model [162,163]. In particular, in ref. [164], the authors consider the three-body decay of the dark vector mesons V into the dark pions π, V →π + + − . In the literature, dark photons are assumed to be heavier than dark mesons in the dark sector (m A /mπ = 3 , m V /mπ = 1.8). The dark photon predominantly produced via Bremsstrahlung in a collision, and their prompt decay produces the dark vector mesons. We map their results into the dark-neutral pion decay. Concerning the production through the off-shell dark photon, we use eq. (4.2) and compute the off-shell dark photon decay rate Γ A through the WZW term [158,159] The production processes of dark-neutral pions through the interactions are depicted in figure 3. Unlike elementary particles, the production of the dark-neutral pions increases following the injected energy due to the derivative couplings. The produced numbers of the dark-neutral pions would be dominated by the UV contribution. The decay rates of the off-shell dark photon (with the fictitious mass m * A ) are given by (4.16) Here, we assume m * A m π for the decay A * → π + π − π 3 . The derivative interaction of dark photons and dark pions lead to the enhancement of the production cross section when a higher energy is injected to the dark sector. The three-pion production is negligible since the production is highly suppressed by high power of m * A /f π and since we use the dark hadronization when the m * A is larger than f π . We discuss the production of the dark-neutral pions via off-shell dark photons at the fixed-target experiments, E137 and SeaQuest. As for E137, the scaling of N A via electron Bremsstrahlung is 2 m −2 A [4]. Since each logarithmic bin of m * A equally contributes to the dark pion production, we cannot identify a certain m * A at which the dark pions are predominantly produced. Hence, we do not factorize the signal number eq. (4.3) as the product of the produced number and the efficiency at a certain m * A . The high m * A contribution to the signal number is suppressed due to the angle acceptance and the loss of the nuclear coherence. Hence, the signal number of dark pions is dominated by the IR contribution.
As for SeaQuest, the dark photons with the mass lighter than 2 GeV are predominantly produced via meson decays and proton Bremsstrahlung. The dark photons via meson decays are less boosted compared to those produced via Bremsstrahlung in a proton 6 Our definition of pion decay constant differs from the literature [159] by a factor 2.

JHEP03(2022)176
beam and have a large angle from the beam axis. The dark photon production via the Bremsstrahlung is dominated near the ρ meson resonance at m A ∼ m ρ , and is nearly constant for m A 500 MeV. It is reasonable to assume that the dark pion production is determined by a virtual dark photon production with a fictitious mass of m * A = m ρ . Since the ρ meson has a broad decay width Γ ρ , we integrate the produced number of A within the range of m ρ ± Γ ρ . Though we compute the signal number of dark pions by use of eq. (4.3) numerically, the π yield is approximated by Here, the first term comes from the left process, while the second term comes from the right process in figure 3. The dark pions can be produced via hadronization since the beam energy of SeaQuest is 120 GeV, which is much higher than the dark dynamical scale. However, the Bremsstrahlung production gets suppressed above 1 GeV due to the proton form factor, and then the Drell-Yan production of dark photons begin to be relevant above 2 GeV. Compared with the Bremsstrahlung at low m A , the A yield via Drell-Yan process is much smaller at high A mass even though the Drell-Yan dominates the production. For the production of the dark pions at the SeaQuest, therefore, we assume that the dark pions are predominantly produced via proton Bremsstrahlung near the ρ-resonance. As for the scaling of N A via proton Bremsstrahlung, we use N A at SeaQuest computed in ref. [21]. Now, let us discuss the sensitivity curves of the visible decay searches for the darkneutral pions. We focus on the process π → A + + − , but the final state A would promptly decay into leptons when the kinetic mixing 10 −4 . In this study, we ignore the decay of final state A to a lepton pair. As in the dark nucleon searches [see eqs. (4.7) and (4.9)], we use approximation formulae for the sensitivity boundaries instead of simulating the production number and the signal efficiency, again.
The first equation gives a lower limit on the boosted decay length of π , while the second gives an upper limit. N 0 LLP and A 0 eff denote the produced number of the long-lived particles and the efficiency in the references, respectively. The produced number N π and the efficiency A eff appears as the ratios to those in the reference model, and hence we only have to care about the parameter dependence of the produced number and the efficiency.
As for the LHC lifetime frontier, we use the same reference to estimate the boosted decay lengths d 0 (p min ) and d 0 (p max ) as in the previous subsection, and use the produced number of long-lived particles N 0 χ in the inelastic DM model for N 0 LLP . We take the momentum p * max to be the same as that of the dark nucleon, which are shown in table 2. Meanwhile, the momentum p * min is different from the dark nucleon case. Since the mass difference of dark-neutral pion and dark photon is larger than that of dark nucleons, the -27 -JHEP03(2022)176 decay products are more energetic. The dark pions leaving signals typically have the momentum of p * geo = m * A /θ det with the fictitious mass m * A which is comparable to the dark dynamical scale and the detector angle θ det with respect to the beam axis. The typical momentum is larger than the threshold momentum as long as the dark dynamical scale is larger than 1 GeV, and then we take p * min = p * geo at MATHUSLA and FASER. As for E137, we take r min = 179 m and r max = 383 m as shown in table 2. We map the result of the visible decay search for dark photon at E137 [4] into the dark-neutral pion search, and N 0 LLP is the number of dark photon produced at the E137. We read out two reference parameter sets ( , m A ) which, respectively, correspond to reference points on the upper and lower boundaries of sensitivity area from ref. [4]. p * max = 20 GeV is given by the maximum momentum of produced dark pions that corresponds to the electron beam energy. Meanwhile, p * min corresponds to the momentum of dark pions that typically comes to the decay volume at E137. The differential production cross section of the dark photon has the collinear singularity, and it is regularized by the dark photon mass [4]. The momentum of the off-shell dark photon peaks near the beam energy, 7 and hence we take p * min = 20 GeV and the production efficiency is not suppressed. We include the angle acceptance at E137 as A eff (p * ) = (1+θ 2 /θ 2 det ) −1 with θ = m * A /p * and θ det 0.004 in the integral in eq. (4.3). Due to the loss of nuclear coherence at the high mass m * A , the contribution to the signal number is suppressed for the high m * A . We do not take into account the decoherence in the integral in eq. (4.3), and hence we take the high-energy cut of the integral to be at m * A 1 GeV. We take the same momentum for the on-shell dark photon, p max = p min = 20 GeV, to evaluate d 0 (p min ) and d 0 (p max ) as for the mapping of the visible decay search for dark photon.
Concerning SeaQuest, we take r min = 5 m and r max = 6 m as shown in table 2. p * max = 120 GeV is the momentum of produced dark pions that corresponds to the proton beam energy. For the SeaQuest sensitivity, we incorporate dark photon production only via proton Bremsstrahlung, and the off-shell dark photons are predominantly produced at the resonance, m * A = m ρ . Meanwhile, the angular scale of the SeaQuest spectrometer is θ det = 0.05 [21]. Hence, the dark pions with the momentum p * min = p * geo = m ρ /θ det 16 GeV predominantly come to the detector at the SeaQuest. We incorporate the geometric acceptance as A geo (p * ) = (1 + θ 2 /θ 2 det ) −1 with θ = m ρ /p * . As for the production through the proton Bremsstrahlung, the differential production cross section peaks at the soft dark photon emission, and hence we take A prod (p * ) = p beam /p * A where p * A is the momentum sum of the final states of Bremsstrahlung production and is constant when the momentum p * is fixed. We map the result of the search for the three-body decays of dark vector mesons at SeaQuest [164] into the dark-neutral pion decay. We take p max = 120 GeV as for the momentum of dark vector mesons. Meanwhile, due to the IR regulator of the soft emission via the proton Bremsstrahlung, we take the typical momentum of dark vector meson leaving signal to be p min = max(m A /θ det , 10 GeV) where 10 GeV is assumed as the minimum momentum in the literature [21]. We read out two reference parameter sets ( , m A ) which, respectively, correspond to reference points on the upper and lower 7 Indeed, the energy distribution of the production cross section implies that the cost for the production of Here, A prod is normalized by the IR regulator p A |max = p beam − m 2 A /p beam .

JHEP03(2022)176
boundaries of sensitivity area from ref. [164], and then we obtain d 0 (p min ) and d 0 (p max ) for the dark vector meson searches at SeaQuest.

Results
As shown in the previous section, there are three kinds of long-lived particles in the composite ADM models: dark nucleons, dark pions, and dark photon. We take the ADM mass to be m N 2 = 8.5/n g GeV. In particular, we take n g = 1 , 2 , 4 , and 8 (m N 2 = 8.5 , 4.3 , 2.1 , and 1.1 GeV). n g = 8 corresponds to the maximum number of generations with which the dark QCD is asymptotically free. The dark QED coupling is assumed to be α = 0.05 , 0.03 , 0.01 , and 7 × 10 −3 , respectively, for n g = 1 , 2 , 4 , and 8 in order to avoid the Landau pole up to the Planck scale. Since we assume m A < m π < 2m A in this article, dark photons are the lightest particle in the dark sector. We place the constraints and the future sensitivities from the dark photon searches in the same plots. First, we discuss searches for the dark nucleon decay, N 1 → N 2 ff , at the LHC lifetime frontier. Figure 4 shows the future sensitivities along with the dark photon searches in the dark photon parameter -m A plane. In this plot, we take the nucleon mass difference to be ∆ N = 0.3 (left) and ∆ N = 0.1 (right). We take different mixing angle θ V in the left and right panels: sin θ V = 5 × 10 −3 (left) and sin θ V = 5 × 10 −2 (right). As for the dark photon constraint, the shaded region is excluded by the existing constraints: the top of figure is already excluded by BaBar, KLOE, and LHCb. The future prospects of the dark photon searches are shown as thin-dashed lines in the figure: Belle-II (brown) and LHCb (orange) on the top of panels.
The thick lines show the sensitivity by the LHC lifetime frontier with the ADM mass to be fixed, MATHUSLA (green) and FASER (blue). We plot the cases with n g = 1 , 2 , 4 , and 8 in each panel, which correspond to the dark nucleon mass to be m N 2 = 8.5 , 4.3 , 2.1 , and 1.1 GeV, respectively. The sensitivity curves for different n g are differentiated from each other by their line types. We use the approximation formulae eqs. (4.7) and (4.9) to draw these curves. Since the dark nucleon mass is fixed, the lifetime of N 1 gets longer as the dark photon mass gets larger; that is accessible at the LHC lifetime frontier gets larger when the dark photon mass increases. The decay mode N 1 → N 2 + A opens when the dark photon mass is lighter than ∆ N m N 2 , and therefore the plots are sharply cut at the left side of the sensitivity region. We assume that the dark photon is the lightest particle in the dark sector, and therefore we cut the sensitivity region at m N 2 = m A .
The mixing angle θ V measures the fraction of darkly charged baryons in the lightest dark nucleon. As discussed in section 3.2, direct detection experiments have put constraint on mixing parameters and θ V since the darkly charged baryons can interact with the SM particles via dark photon. The direct detection bounds are placed the left-top region of each panel, and the region above the lines are excluded. The four line types correspond to that of sensitivity plots: the mass of ADM is 8.5   [3,8], KLOE [9,[80][81][82], and LHCb [126,165]; (left-bottom region) the fixed-target experiments νCal [144,145] and CHARM [141,142]. The future sensitivity of Belle-II (LHCb) to visible dark photon decay is shown as brown (orange) dashed lines on middle parts of the panels [25,27]. The gray diagonal lines shows the direct detection bound on the composite ADM, and the parameter space above the line is already excluded. The sensitivity of MATHUSLA and FASER sharply cut at the left side of sensitivity because the decay into the on-shell dark photon opens.
The produced particles leaving signals are less boosted for the MATHUSLA sensitivities than for the FASER sensitivities since MATHUSLA will be located in the off-axis direction from the LHC beam line. Thus, MATHUSLA has the sensitivity to the longer lifetime, namely smaller , than FASER even though the both detectors will be located at the similar distance from the LHC collision point. Since the DM is produced via the Drell-Yan mechanism, the produced number enhances for the lighter DM. Besides, the DM is produced more preferably in the forward direction to the beam line as the DM gets lighter, and hence the FASER sensitivity area gets much wider than the MATHUSLA one. The dark nucleon searches at the LHC frontier have the sensitivity to 10 −4 . In other words, some portion of the sensitivity to the dark nucleon visible decay is comparable to the future sensitivity at the LHC-b and Belle-II that searches for prompt decay of dark photons. In the plots, the upper and lower sensitivity curves sharply cross each other at the left-bottom point. This is because we use the approximation formulae the sensitivity curves, eqs. (4.7) and (4.9), and then it is expected that the shape of sensitivity curves of MATHUSLA and FASER will be round on the edge. Most of the parameter space where the approximation formulae are robust has been excluded by the current bound on the visible dark photon decay. The dedicated analysis may lower the lower sensitivity curves since our approximation formula provides the conservative lower bound for multi-GeV dark nucleon (see appendix C).
We show the parameter dependence in the figures; in particular, we change the mixing angle θ V , the mass difference ∆ N , and the mass of DM m N 2 . In the following, we discuss how the sensitivity curves change when we make them larger. The lifetime gets shorter as the values of these parameters get larger. Since the upper boundary of the sensitivity area is mainly determined by the lifetime, the smaller will be explored in order to fix the lifetime. On the other hand, the lower boundary of the sensitivity area is basically determined by the product of lifetime and production cross section. The production cross section does not significantly depend on θ V and ∆ N , but does depend on 2 . The lower boundary also gets lower in as the values of these parameters get larger, but not so much as the upper boundary. As a result, the sensitivity range of shrinks for large values of parameters. For instance, figure 4 shows the sensitivity area gets shrinking for the heavier DM mass (smaller n g ). In particular, sensitivities disappear for n g = 1 except for the right top panel of figure 4. Once ∆ N m N 2 exceeds the dark photon mass, the decay mode with the on-shell dark photon opens. The sensitivity range of m A also shrinks for large ∆ N . We also show different choices of parameters in figure 11 in section D.
We take the multiplicity of dark nucleon production to be a similar to that of nucleons in the SM near J/ψ threshold. In the SM, pions are much lighter than nucleons, and hence it would be expected that pions are likely to be produced more than nucleons. In our model, on the other hand, the dark pion mass can get closer to the dark nucleon mass compared to the SM case. The multiplicity of dark nucleons can be similar to the dark pions, and then the sensitivity range would be enhanced. We discuss the change of multiplicity in section D.
We note that the fixed-target experiments, such as E137 and SeaQuest, would also be available to searching for the dark nucleon below a few GeV (cf., sub-GeV inelastic DM  searches [21,146]). We do not investigate this possibility for the lighter dark nucleons further in this study. Figure 5 shows the future sensitivity and the existing constraint on the dark pion three-body decay, π 3 → A + ff . There are three dimensionful parameters associated with -32 -

JHEP03(2022)176
the dark pion decay: dark photon mass m A , dark pion mass m π , and dark pion decay constant f π . We assume the dark pion mass in m A < m π < 2m A in order that the three-body decay dominates the dark pion lifetime, in particular we take m π /m A = 1.9 (solid lines) and 1.3 (dashed lines). Since the dark dynamical scale is fixed by the ADM mass, the decay constant can be scaled as a function of dark nucleon mass, The shaded region on the top-left of each panel is excluded by the E137 bound: we take m π /m A = 1.3 (with dashed-line boundary) and m π /m A = 1.9 (with solid-line boundary). E137 uses the electron beam of 20 GeV, and hence it is hard to produce heavy particles of multi-GeV through the electron Bremsstrahlung at E137. We take the fictitious dark photon mass to be less than 1 GeV where the form factor of target nucleus will be O(1). When the energy injection to the dark sector is above the dark dynamical scale, it is possible to produce dark pions through dark hadronization. We use the hadronization for the dark pion production eq. (4.14) with f π = 0.1 GeV, while use the production via the off-shell dark photon eq. (4.15) for other f π . The dark hadronization must produce six dark pions at one hadronic event due to the isospin symmetry and the multiplicity of the dark pion production n π = 2. The exclusion area for f π = 0.1 GeV is sharply cut at 0.2 m A /m π GeV where the total mass of the final states via the dark hadronization, 6m π , equals to the high mass cut 1 GeV as discussed in section 4.2.
The magenta lines show the future sensitivity of SeaQuest to the three-body decay, while the cyan and the yellow lines show the future sensitivity at the LHC lifetime frontier, FASER and MATHUSLA, respectively. At SeaQuest, we focus only on off-shell dark photon production via Bremsstrahlung in a proton-nucleus collision though it would be also important to take into account the dark photon production via Drell-Yan processes above a few GeV. Since the off-shell dark photon production is sharply dropped above the fictitious mass m * A m ρ , which we use to estimate the signal number of dark pions in eq. (4.3), the SeaQuest abruptly loose their sensitivity to the dark pion decay near a few hundred MeV. Similarly to the E137 constraint, we use the hadronization for the dark-neutral pion production eq. (4.14) with f π = 0.1 GeV, while use the production via the off-shell dark photon eq. (4.15) for others.
As shown in figure 5, the sensitivity curves of SeaQuest and the existing constraints by E137 drastically change since the production channel of dark-neutral pions via dark hadronization opens as the injected energy is higher than the dark dynamical scale. Even in the SM, the production cross section with the neutral pion final states drastically change near the dynamical scale. As for the charged pions in the SM, the production cross section -33 -

JHEP03(2022)176
gradually increases due to the ρ-meson broad width as the injected energy to the hadronic sector increases. Meanwhile, the production of neutral pions opens at the ω-meson threshold [89,160,161]. Due to the narrow width of the ω-meson, the production cross section sharply increases at the dynamical scale. We do not go further into the hadronization in this study since we are agnostic about the hadronization.
At LHC lifetime frontier, the dark pions are produced via three dark pion production with an injected energy to the dark sector above the dark dynamical scale. It is required to inject more energy in order to produce dark pions when 6m π exceeds the dark dynamical scale. In this study, we assume 6m π ≤ m N 2 instead of dedicated analysis of hadronization in the dark sector, and hence the sensitivity region is sharply cut at the threshold As we discuss in the dark nucleon searches, MATHUSLA has the sensitivity to the longer lifetime, namely smaller , than FASER. Similarly to the dark nucleon searches, the dark pion searches at SeaQuest and LHC lifetime frontier are comparable to the future sensitivity of prompt decay searches for dark photons.  ). Throughout this study, we take the dark QED coupling not to be diverged up to the Planck scale. We can take a larger coupling when some of generations are decoupled from the low-energy theory and the dark QED is unified into Yang-Mills theory with a larger gauge group which is asymptotically free at the high-energy scale. 8 In that case, both of the decay rate and the production cross section of the dark-sector particles increase, and hence we can explore more smaller kinetic mixing .
Furthermore, the dark hadrons produced at the visible decay searches are assumed to be made only of one of the n g generations in this study. In other words, we assume that all dark quarks except for the lightest generation have masses larger than O(10) GeV or more. The production rate of dark hadrons enhances due to transition via the dark strong force from other dark hadrons when all n g generations are produced. In this case, the production rate is naïvely expected to be multiplied by n g , and hence it makes the lower boundary of the sensitivity curves a factor of n 1/4 g smaller than what we have shown.

Conclusion
Composite ADM models with dark photon have multiple particles of GeV or sub-GeV scale: dark photons, dark pions, and dark nucleons. Each of particles plays significant roles in the composite ADM framework: the lightest dark nucleon is the very ADM, the strong 8 The dark nucleon mass is determined by the degree of freedom in the dark sector at the decoupling temperature of the intermediate-scale portal interaction eq. (2.2). Therefore, the dark quarks should be decoupled (just) below the decoupling temperature scale. Meanwhile, in order not to change the degrees of freedom in the dark sector, the dark QED is unified into Yang-Mills theory above the decoupling temperature scale. annihilation into dark pions depletes the symmetric component of dark nucleons, and the dark photons release a huge entropy in the dark sector to the visible sector. The dark sector connects to the visible sector via the dark photon portal only with small kinetic mixing at a low energy, and thus the dark particles tend to be long-lived and leave visible signals. The lightest dark nucleon must have the mass of GeV to have the correct relic abundance in the composite ADM framework, and therefore it is certainly a good target of the LHC lifetime frontier.
In this study, we have focused on the case that the lightest dark nucleon mainly consists of dark neutron but is slightly mixed with dark proton due to U(1) D breaking. The constraint from the direct detection on the kinetic mixing gets milder thanks to the nucleon mixing. In this case, the heavier nucleon that is mostly composed of dark proton can decay into the lightest dark nucleon through three-body decay. Figure 4 summarizes the dark nucleon searches at the LHC lifetime frontier, FASER and MATHUSLA. They will explore dark photon mass above sub-GeV and kinetic mixing 10 −4 -10 −3 for the DM of multi-GeV.
Meanwhile, dark pions are lighter than dark nucleons, which means that the dark pions have the mass of sub-GeV in the composite ADM framework. In this study, we have considered dark pions with the mass in m A < m π < 2m A that is optimized for the visible signals from dark pions via three-body decay π 3 → A + ff . E137 puts the most strong constraint on the decay, and SeaQuest have a great sensitivity to the decay in near future exploring dark photon mass below GeV and kinetic mixing 10 −5 -10 −4 Besides, we have discussed the dark pion searches at the LHC lifetime frontier. The dark pions have been assumed to be produced via hadronization in the dark sector at the LHC. Similarly to the sensitivity of SeaQuest, we have found that the LHC lifetime frontier has the potential to explore the kinetic mixing of 10 −4 . Dark strong dynamics naturally provides rich structure in the dark sector. Hence, we have various decay signals from dark hadrons unlike, for example, the visible decay of dark photon. Interestingly, the visible signals from dark hadrons and dark photon will be tested at various experiments in near future: dark nucleons at LHC, dark pions at LHC and SeaQuest (fixed-target experiments), and dark photons at Belle-II and LHCb.
In this study, we use the approximation formulae to estimate the sensitivity boundaries. In order to precisely determine the future sensitivity, we need to understand the produced numbers of dark hadrons and the efficiencies. In particular, dark hadrons are assumed to be produced via the off-shell dark photon or via the dark hadronization in this study. The dark-sector particle production through the off-shell A has been less studied compared to the production through the decay of on-shell A . We evaluate the produced numbers of dark hadrons via the off-shell dark photons by use of the on-shell dark photon production at the fictitious mass m * A and by considering the dominant contribution to the numbers. It is important to study the off-shell production in more detail, in particular by simulating the produced numbers of dark hadrons at LHC and at the fixed-target experiments, for the future sensitivity. Besides the dedicated analysis for the off-shell production, it is also required to understand hadron physics with different parameters from the SM ones. We have naïvely used empirical values of hadron production in the -36 -

JHEP03(2022)176
SM in order to estimate produced numbers of dark hadrons. It is worth studying more dedicated analysis with inclusion of dark hadron physics (e.g., R-ratio of dark hadrons, multiplicity of dark hadrons, etc.), but we leave them for future study.
In this study, we have focused on a specific mass spectrum of the dark pion, m A < m π < 2m A , which is relevant for visible decay searches. For the dark pions with 2m A < m π , produced dark pions promptly decay into dark photons through the anomaly-induced interaction, π 3 → A A , and then we lose the future sensitivity and the existing constraints from fixed-target experiments. When the dark pions are the lightest particle in the dark sector, m π < m A , dark pions can decay into the visible particle with lifetime longer than 1 s. Instead, dark photon would dominantly decay into dark pions providing invisible signals of dark photons at the fixed-target experiments such as LDMX [22,23]. Furthermore, since the dark sector consists of the confining gauge dynamics, dark vector mesons would exist with the mass of dark dynamical scale. If the dark vector mesons have a broad resonance similar to the SM ρ mesons, the production of darkly charged pions would change due to the broad resonance and kinetic mixing between dark photons and dark ρ mesons. We leave the long-lived particle searches with different mass spectra for future work.

A Linear sigma model
We analyze chiral symmetry breaking and scalar dark mesons by use of the linear sigma model (LSM) that is a low-energy effective theory of two-flavor dark QCD. The Lagrangian density of the LSM is constructed to be invariant under the flavor transformation, SU(2) L × SU(2) R ×U(1) V , which reflects the symmetries of the dark QCD. U(1) A is explicitly broken by an axial anomaly term. Let Φ be an LSM field that is an 2 × 2 matrix scalar field. Φ is identified with the dark-quark bilinear as Φ ij ∼q j q i with the flavor subscripts i, j. According to this identification, Φ transforms under the flavor transformation, SU(2) L × SU(2) R × U(1) A , as follows. determine C and P transformation of Φ as follows.
Here, we define the discrete transformation for the field in the charge basis, which is denoted by the hatted fields. We decompose the matrix field Φ into scalar fields σ a and pseudo-scalar fields π a as where T a (a = 0, 1, 2, 3) are a unit matrix and SU(2) generators that satisfy f abc is totally antisymmetric in all indices, while d abc is totally symmetric. They are just structure constants of SU(2) when a, b, c = 0, while we have f ab0 = 0 , d ab0 = δ ab . The Lagrangian density is In this sub section, we do not include the gauge interactions of the LSM field. Here, the second term is the potential for Φ that is given by We assume parameters λ 1 , λ 2 , and µ 2 to be positive. The last term breaks U(1) A , and the coefficient c is a complex number with mass dimension two in general. Throughout this paper, we fix the coefficient c to be positive by U(1) A rotation of Φ.
The third term of the Lagrangian reflects the explicit breaking of the chiral symmetry by the dark quark masses. We refer to H ≡ j a T a as a source matrix that is proportional to the mass matrix of the current dark quarks. It would be sufficient to consider the linear term in H as far as the dark quark masses are sufficiently small and its perturbation works well.
The Lagrangian is rewritten in terms of the components scalars, σ a and π a , as follows.

JHEP03(2022)176
σ 0 and π a (a = 1, 2, 3) may obtain the dominant vacuum expectation value (VEV) since they apparently have negative mass squared. We simply set π a that obtains its VEV to be π 3 by use of the vectorial flavor rotation SU(2) V of Φ. We parametrize masses of scalars after the symmetry breaking as follows.
Here, letters with bars denote their VEVs, and the higher terms are abbreviated. Let us compute the mass spectra under some assumptions. We can generally diagonalize the source term by use of SU(2) L × SU(2) R transformation of Φ.
We can utilize U(1) A transformation to remove the overall phase α H , but the phase of the anomaly term c appears, again. To avoid the phase of c, we keep α H unless we assume the source term to be real. The kinetic term and the potential V (Φ) is invariant under Φ → −Φ as the number of flavors N f is even. Therefore, we can set j 0 > 0 without loss of generality for N f = 2.

A.1 Decay constant
Let us discuss the decay constant in the LSM. In particular, we clarify the normalization of the decay constant by use of the commutation relations of the charges and partiallyconserved axial current (PCAC) relation. The infinitesimal axial-vector (L = R † = e iα a T a ) and vector (L = R = e iβ a T a ) transformations define their Noether currents, J µ Aa and J µ V a .
Here, we define Φ ≡ φ a T a . The PCAC relation gives where f ab is the decay constant. In terms of the component fields, σ a and π a , the axial current is J µ . When the sigma fields get their VEVs, σ a =σ a , the PCAC relation gives the decay constant as follows.
In particular, when the VEV is universal,σ a =σ, the decay constant is diagonalized: We can also define the currents in terms of quarks denoted by Q: This quark axial current annihilates the one-pion state, and it also defines the pion decay constant.

JHEP03(2022)176
In the SM, the observation of the charged pion decay, π + → µ + ν, determines the decay constant. In this normalization, we have f π = 92 MeV from the SM pion decay. This normalization provides a commutation relation of the charges as follows.
In the LSM, the axial charge and the vector charge are defined by These charges obey the following commutation relation.
The charges in the LSM satisfy the same algebra as the charges defined in terms of quark currents. Thus, we conclude that f ab defined in eq. (A.13) follows the same normalization as the decay constant defined by quark currents.

A.2 Universal source
In this subsection, we compute the spectrum with the source term H = e iα H j 0 T 0 (j 0 > 0). First, we consider the case with the real source (α H = 0). We set σ a =σ 0 δ a,0 and π a = 0. The tree-level effective potential is given by Since the potential is tilted in the presence of j 0 > 0, the stationary point condition for this potential gives the VEV up to O(j 0 ) as .
We define the decay constant f ≡σ 0 . We get the mass spectrum up to O(j 0 ) as follows.
Here, the subscript i runs 1 , 2 , and 3. There is no off diagonal components of mass matrices, and thus the mass spectrum is already diagonalized. Since the source j 0 is proportional to the current quark mass, π i (i = 1, 2, 3) are identified to be the pseudo Nambu-Goldstone bosons. We can identify σ 0 to be the so-called σ meson, and π 0 to be the so-called η meson. The U(1) A anomaly provides the sizable contribution to the mass of π 0 , and therefore π 0 is heavier than the other pseudo scalars.

JHEP03(2022)176
When the phase α H is turned on, the potential is also tilted to the π 0 direction. Therefore, π 0 acquires its VEV that proportional to the source, and we set σ a =σ 0 δ a,0 and π a =π 0 δ a,0 . The effective potential is given by As we mentioned above, we can restrict the parameter space to j 0 cos α H > 0 by flipping the sign of Φ without loss of generality. The stationary point condition for the potential gives their VEVs as follows.
Since both σ 0 and π 0 get their VEVs, there is mass mixing among σ a and π a at the order of j 0 . This leads corrections of O[(j 0 ) 2 ] to the mass eigenvalues. At the leading order of j 0 , we have the similar mass spectrum of scalars to the one without the phase α H .
We define the decay constant f ≡σ 0 , again.

A.3 Isospin-violating source
In this subsection, we compute the spectrum with the source term H = e iα H j 3 T 3 (j 3 > 0). As with the previous subsection, we first consider the case with the real source (α H = 0). Since the potential is tilted in the σ 3 direction, σ 3 must have its VEV that is proportional to the source. However, the dominant VEV is determined irrespective to the source term since the mass squared of σ 3 is positive. As discussed in the introduction of this appendix, there are two possible dominant VEVs: one isσ 0 and another isπ 3 . First, we discuss the former case; σ a =σ 0 δ a,0 +σ 3 δ a, 3 and π a = 0. The effective potential for this vacuum choice is .

JHEP03(2022)176
We note thatσ 3 is positive for j 3 > 0 and is proportional to j 3 . There is no linear correction of j 3 to the scalar masses since the source j 3 appears inσ 0 from the second orders of j 3 . The pseudo Nambu-Goldstone bosons are massless up to O(j 3 ). The scalar mass spectrum Since there is the mass mixing of O(j 3 ) between π 0 and π 3 , we have to take care of corrections of O[(j 3 ) 2 ]. In particular, since the dark pions are massless even at the order of O(j 3 ), it is not obvious if the potential for pseudoscalars is stabilized. The lightest pseudoscalar mass is indeed tachyonic: and hence, the lightest pseudoscalar developes its VEV. We conclude that the choice of VEVs, σ a =σ 0 δ a,0 +σ 3 δ a, 3 and π a = 0, does not provide a proper vacuum. Instead, we discuss the latter case; σ a =σ 3 δ a,3 and π a =π 3 δ a, 3 . The effective potential for this vacuum choice is given by From the stationary point condition of this potential, we determine the VEVs at the order There is no linear correction of j 3 to the scalar masses since the source j 3 appears inπ 3 from the second orders of j 3 . The pseudo Nambu-Goldstone bosons are massless up to O(j 3 ). The scalar mass spectrum up to O(j 3 ) is given by Here, we take f =π 3 as with the previous case: This spectrum coincides with that with the previous vacuum choice up to O(j 3 ). As with the previous case, we have to care about correction up to O[(j 3 ) 2 ] since we have massless scalars. There is mass mixing between σ and π since both σ and π get their VEVs in this vacuum choice. Every massless scalars obtains the mass of O[(j 3 ) 2 ]: Here, j 3 andσ 3 have the same sign, and therefore all massless scalars get positive mass of O[(j 3 ) 2 ]. We conclude that the vacuum choice, σ a =σ 3 δ a,3 and π a =π 3 δ a,3 , provides a proper vacuum up to O[(j 3 ) 2 ] for the isospin-breaking source H = j 3 T 3 (j 3 > 0).

JHEP03(2022)176
We note that mass spectrum eqs. (A.35) and (A.36) is similar to the spectrum eqs. (A.26) and (A.27) with α H = π/2. In particular, the spectrums coincide with each other by replacing π 3 ↔ σ 0 and π 0 ↔ σ 3 . The universal mass with α H = π/2 implies the source term is pure imaginary. Since we can rotate the LSM field by a specific flavor rotation: mass spectrums have correspondence between real and pure imaginary sources. This rotation implies that the mass spectrum of pure imaginary source j 3 coincides with eqs. (A. 22) and (A.23). Let us confirm this correspondence explicitly by compute the spectrum with the complex isospin breaking source H = e iα H j 3 T 3 (j 3 > 0). We again take the same vacuum σ a =σ 3 δ a,3 and π a =π 3 δ a,3 as before.
This determines the VEVs from the stationary point condition of the potential as follows: These VEVs give the mass spectrum of scalars up to O(j 3 ) as follows: When we take α H = π/2 (pure imaginary j 3 source), the mass spectrum coincides with the mass spectrum with the real j 0 source (see eqs. (A.22) and (A.23)).

A.4 Mixture
In the last of this section, we consider more generic source term. In general, quark masses violate the isospin symmetry, and hence the source term in the effective theory should reflect the violation. In particular, the isospin violation in the mass spectrum is important for the pion mass difference in our study. The diagonalized quark mass matrix corresponds to the following source term in the LSM: In the presence of two sources j 0 and j 3 , the vectorial flavor symmetry SU(2) V is also broken. The vacuum choice should reflect the breaking, and hence we assume the vacuum as follows: σ a =σ 0 δ a,0 +σ 3 δ a,3 , π a =π 0 δ a,0 +π 3 δ a, 3 .

JHEP03(2022)176
In the following, we find analytic formulae of the mass spectra. We compute the spectrum by treating either j 0 or j 3 as a dominant source and by taking into account another as perturbation, and then we check the absence of tachyonic modes. First, we consider the case with j 0 j 3 . From the dominant source j 0 , the VEVσ 0 is given by eq. (A. 21), and mass spectrum is given in eqs. (A.22) and (A.23). The potential for σ 3 (mass: 2c+λ 2 f 2 with f 2 = λ −1 (c+µ 2 )) is tilted by the source term j 3 , and therefore σ 3 gets its VEV asσ 3 j 3 (2c + λ 2 f 2 ) −1 as far as the potential is mainly lifted up by the mass term. Thus, the VEVs in the case with j 0 j 3 are given bȳ The presence ofσ 3 leads to the modification of the mass spectrum of the order O[(j 3 ) 2 ], in particular the mass difference between π 3 and π 1,2 . π 0 -π 3 mixing arises fromσ 0 andσ 3 , and then the mass eigenvalues are . (A.47) Here,π 0 andπ 3 are mass eigenvalues of the π 0 -π 3 system. The presence of j 3 leads to the mass difference among π a (a = 1, 2, 3) though they form the isospin triplet. We note that the isospin violation in the pion mass arises of order The condition that the lightest pion has a positive mass squared puts the upper bound on j 3 : j 3 max ≥ j 3 where Next, we consider the case with j 3 j 0 , i.e., we treat j 3 as the dominant source. Eq. (A.34) gives the VEVs when the source term contains only j 3 . Both σ 3 and π 3 get their VEVs in this case since there is the mass mixing among σ and π. Therefore, both σ 0 and π 0 get their VEVs via j 0 though j 0 makes the potential only for σ 0 tilted. The mass matrix for σ 0 and π 0 has components up to O[(j 3 ) 2 ] as where we define f 2 = λ −1 (c + µ 2 ), again. The mass basis and the mass eigenvalues are defined by where the mixing matrix R and the mass eigenvalues m 2 ± are

JHEP03(2022)176
Both s − and s + get their VEVs, and hence we obtain the VEVs in the original basis as follows.
(A. 53) In the limit of small sources, we obtain the VEVs:

A.5 Higgs-scalar meson mixing
In this subsection, we consider back-reaction from the VEV of the dark Higgs φ D . The potential including the LSM field and φ D is given by where each part of the potential is where j a T a is a diagonalized source term that originates from the dark-quark mass matrix. V mix denotes the mixing term of the LSM field and φ D Here, v is the VEV of the dark Higgs that is determined only by the Higgs potential, v = 2µ 2 φ /λ φ , and j Y and Φ D are the matrix-form coupling and the dark Higgs field defined by j Y originates from the Yukawa coupling between the dark quarks and the dark Higgs. j 1 and j 2 are, respectively, proportional to Yukawa couplings y 1 and y 2 in eq. (2.5), and j a (a = 1, 2) have mass dimension of three. A combination Φ D j Y has only off-diagonal entries, and gives a source term of the off-diagonal component of Φ after φ D gets its VEV. We define the Higgs-dependent source term H as follows: where

JHEP03(2022)176
When the source term is hermitian, H = H † (i.e., j 1 = j 2 ), the source matrix H is factored out trH(Φ + Φ † ). In this case, the dark Higgs potential is tilted if the charged scalar meson obtains its VEV in the charged basis. The Lagrangian possesses the parity invariance, Φ → Φ † . As shown in previous subsection, the dominant VEVs of the LSM field areσ 0 and π 3 in mass basis when parity is conserved. When the system possesses the SU(2) V flavor symmetry, no U(1) D -charged sigma meson gets its VEV in charge basis fromσ 0 . Due to the parity invariance, again, no U(1) D -charged sigma meson obtains its VEV in charge basis fromπ 3 in mass basis. Therefore, we conclude that we have no back-reaction from the LSM field to the dark Higgs VEV when we have P -invariance and SU(2) V global flavor symmetry. We therefore consider the following cases that the source term is not hermitian.
Since the quark mass matrix is diagonalized before the U(1) D breaking, the original source term j a T a is diagonalized and we parametrize the source term as follows.
We consider the first case that H = J 2 U 2 A with an SU(2) matrix U A and a real-positive parameter J > 0. The matrix U A is parametrize phase α A and angle θ A as follows: By comparing with eq. (A.63), we determine parameters as follows.
The parameters depend on the Higgs field as follows.
we can diagonalize the source term. Φ ≡ U A ΦU A refers the LSM field in this basis. In this basis, the source term is proportional to unity, and therefore Φ obtains its VEV proportional to unity: Φ =σ 0 1, and hence the VEV of Φ has an off-diagonal entry that provides U(1) D violation.

JHEP03(2022)176
At the LSM potential minimum, we have a linear term of J and this gives a back-reaction potential to the dark Higgs potential.
We define the back-reaction potential as V BR (φ D ), and expand φ D around the VEV v, where f = λ −1 (c + µ 2 ) denotes the decay constant. Once we assume that the shift of the VEV is negligible, the pion mass is given by the source J with φ D = v, m 2 π J(v)/f . This term tilts the Higgs potential, and then the true Higgs VEV is shifted from v, v D = v + δv. We determine δv in a perturbative manner of the source. For simplicity, we take the sizes of j 1 and j 2 are the same, but the signs are opposite to satisfy eq. (A.66): and then the shift of the Higgs VEV is given by We use the relation between the source and the pion mass, m 2 π j 2 0 + j 2 /f √ 2j/f , and use the Gell-Mann-Oaks-Renner (GOR) relation, m 2 π f 2 m q qq in the second equation. When j 0 j, the potential is and then the shift of the Higgs VEV is given by In the second equation, we use m 2 and then the shift of the Higgs VEV is given by We use m 2 π j 2 0 + j 2 /f j/f and the GOR relation, again.
-47 - (2) matrices U L and U R , and a real-positive parameter J > 0,
The parameters depend on the Higgs field as follows.
By use of the SU(2) L × SU(2) R rotation, U L and U R , we can diagonalize the source term. Φ ≡ U L ΦU † R refers the LSM field in this basis. In this basis, the source term is proportional to T 3 , and therefore Φ obtains its VEV proportional to T 3 : as shown in appendix A.3, both σ 3 and π 3 get their VEVs, Φ = (σ 3 + iπ 3 )T 3 , and hence the VEV of Φ has an off-diagonal entry that provides U(1) D violation.
Since the VEVπ 3 = f does not contain the linear term of J, the pion mass is proportional to J 2 , where f = λ −1 (c + µ 2 ) denotes the decay constant, again. We assume that the shift of the VEV is negligible, and then the pion mass is given by the source J with φ D = v. At the LSM potential minimum, we have a quadratic term of J and this gives a back-reaction potential to the dark Higgs potential.
Similar to the previous case, we define the back-reaction potential as V BR (φ D ), and expand This gives a shift of the Higgs field as follows up to the second order of sources.
We consider more generic source term that is characterized by T 0 and T 3 . This is a combination of the aforementioned two cases. We consider (2) matrices U L and U R , which we use in the previous case, and a real-positive parameter J 0 , J 3 > 0. We take α L = α R ≡ α V for simplicity, again.
The parameters depend on the Higgs field as follows.
By use of the SU(2) L × SU(2) R rotation, U L and U R , we can diagonalize the source term. Φ ≡ U L ΦU † R refers the LSM field in this basis. As with the previous case, while the VEVs are diagonal in the mass basis, there are off-diagonal entries in Φ .
It is challenging to find an analytic formula that is valid whole parameter range. In appendix A.4, we find the approximate solutions where each of sources is larger than another and we treat the smaller one as the perturbation. When J 0 dominates the source term, the VEVs are given by eqs. (A.45) and (A. 46), At the LSM potential minimum, we have a linear term of J 0 and this gives a back-reaction potential to the dark Higgs potential.
This back-reaction potential gives the same results with the case of H = J 2 U 2 A . On the other hand, when J 3 dominates the source term, the VEVs are given by eqs. (A.54) and (A.55),

JHEP03(2022)176
As with the case of H = JU † R T 3 U L , at the LSM potential minimum, we only have a quadratic term of J 3 and this gives a back-reaction potential to the dark Higgs potential.
This back-reaction potential gives the same results with the case of H = JR † T 3 L.

B Hadron spectrum
In the previous appendix, we discuss the generic mass spectrum in the LSM without corrections from QED. In this appendix, we incorporate the QED correction into the mass spectrum. While the QED corrections usually arise as the loop corrections, we treat them as effective interactions by regarding the charge matrices as spurion fields.

B.1 Pion mass matrix
Let us consider meson spectrum that originates from the universal source discussed in appendix A.2. The meson masses are Here, we assume that the source j 0 is positive. When the quark mass basis differs from its charge basis, we obtain the charge basis by the flavor rotation. It is expected that the charge basis of mesons is identified with the mass basis rotated by the same flavor rotation. Let Φ andΦ be the LSM fields in the mass basis and in the charge basis, respectively.
where the flavor rotation matrices U L and U R that relate the quark mass basis to its charge bases. In particular, when we take U L = U R = U V and their phases to be α L = α R = 0, the charge neutral pionπ is a linear combination of π 1 and π 3 , Since the charge neutral pion has a coupling to dark photons via chiral anomaly, pions can decay through the (off-shell) dark photons through the mixing. Once we have U(1) D coupling, we can also add the charge masses for the dark pions. In general, the masses in the mass basis are given by

JHEP03(2022)176
Here, q L and q R are the charge matrices that are defined in eq. (2.14), and Π is the matrixform pion fields in the mass basis.
Compared to the chiral perturbation theory, Π † Π is not constant in the LSM, and then we can add the second and the third terms. Indeed, these terms are important for giving the correct charge masses. When L and R that relate the quark mass basis to its charge bases are unity, the mass basis coincides with the charge basis. In this case, the charge matrices are the same q L = q R =q = diag( 2 3 , − 1 3 ) and π 1 and π 2 are the charged pions, then we obtain In the last line, we determine the coefficients in order to correctly give charge masses to the charged pions, c 2 = c 3 = −2c 1 ≡ c π Q /2. Next, we take into account isospin breaking terms. For simplicity, we take the flavor rotation matrices U L = U R ≡ U V , and then the charge matrices are the same In this case, the mass terms of pions are given by where we define m 2 D = 4πα c π Q f 2 π . m 2 π 3 and m 2 π 1,2 are defined in eq. (A.47), which include the isospin violation. The mass eigenvalues are where c α = cos α and s α = sin α, and the mixing angle β V is defined by When there is no U(1) D correction, m D = 0, β V is zero. If there is no isospin violation in mass basis in advance, m 2 π 3 = m 2 π 1 , the degeneracy of the dark pions is resolved only by -51 -

JHEP03(2022)176
the U(1) D correction. This means that the charge basis and the mass basis coincide with each other, and indeed β V = −2θ V as m 2 π 3 = m 2 π 1 . The mass eigenvalues are In this mass basis, U(1) D gauge interaction of pions is given by

B.2 Nucleon masses
We now consider the mass spectrum of the light dark nucleons. The dominant nucleon masses arise from the interaction with the LSM fields, while the global symmetry allows the vector-like mass term that originates from the quark mass term.
Here, we construct the hadron effective theory after we rotate the quark basis into their mass basis, and hence the vector-like mass term is proportional to the diagonalized quark mass matrix. Even though the vector-like mass term is not dominant contribution, it is important for the nucleon mass difference. The large-N C scaling of a coupling g is g = g πN N ∼ 4π √ N C . We assume that the LSM field obtains its VEV to be diagonalized, Φ = f T 0 + δT 3 . The nucleon masses are given by eq. (2.24). We can also include U(1) D correction to the nucleon mass, 14) The charge matrices are defined in the text, eq. (2.27). We discuss the nucleon mass difference in the SM, proton-neutron mass difference, by use of the LSM. In the SM, N 1 and N 2 correspond to neutron and proton, respectively. First, we consider the mass difference originates from isospin-violation by quark masses. In the SM, the mass difference of quarks is smaller than their sum, hence we use VEVs given by eq. (A. 45).
In the second equality, we rewrite the mass difference in terms of quark masses and hadron masses by use of the following relations.
Since the source terms originate from the quark mass terms, the ratio j 3 /j 0 is determined by the quark mass difference and their sum. The second relation is obtained from the sum Our fitting with the results of ref. [146] for MATHUSLA (left) and FASER (right). We take the parameters to be m A = 3m χ and α = 0.1, which are the same as those in ref. [146]. The mass difference parameter, ∆, is taken to be ∆ = 0.1 (solid), ∆ = 0.05 (dashed), and ∆ = 0.03 (dotted). The green (blue) lines show the future sensitivities at MATHUSLA (FASER) taken from ref. [146]. The light green (cyan) points in the plot are the reference points of MATHUSLA (FASER) that we have used in the text. The black lines show our approximate sensitivity given by eq. (4.7) and eq. (4.9). of nucleon masses as quark masses are negligible. Numerically, we take the SM values of parameters [90]:

C Fitting
In this appendix, we demonstrate the validity of our approximation formulae eq. (4.7) and eq. (4.9) by applying the formulae to the specific models. As we discussed in the text, as for the LHC lifetime frontier, FASER and MATHUSLA, we apply our formulae to the inelastic DM model [146]. On the other hand, as for the fixed-target experiments, SeaQuest and E137, we apply them to the dark vector meson search [164] for SeaQuest, and to the dark photon search [4].

C.1 MATHUSLA and FASER
We show our fitting in figure 7 by use of our approximation formulae, eq. (4.7) and eq. (4.9). The MATHUSLA sensitivity plots are shown as green lines in the left panel, while the -53 -

JHEP03(2022)176
FASER plots are shown as blue lines in the right panel, which are read from ref. [146]. We take the dark photon mass to be m A = 3m χ with the dark matter mass m χ and the U(1) D coupling to be α = 0.1 that are same as the parameters taken in the literature. The on-shell A is produced via the SM meson decay, the Bremsstrahlung, and Drell-Yan at the LHC, and its decay A → χ 1 χ 2 produces the long-lived particle χ 2 . The mass difference of the inelastic dark matter, defined by ∆ ≡ (m χ 2 − m χ )/m χ with the masses of the excited state m χ 2 , is taken to be ∆ = 0.1 (solid), ∆ = 0.05 (dashed), and ∆ = 0.03 (dotted). We assume that the long-lived particles are produced only via Drell-Yan production at the LHC. We do not take into account Z-resonance for the production, and hence the peaky behavior around m χ 30 GeV (corresponding to m A m Z ) are not reproduced in these figures. We take into account the difference of m A -dependence of N A of χ 2 , which are scaled as 2 m −4 A in forward direction and 2 m −2 A in total, as we mentioned in the text. We take into account the efficiency due to the energy threshold as discussed in the text, A eff = A geo A prod . p min is the larger momentum between p geo and p thr ; p thr = 1 TeV(0.1/∆) for FASER and p thr = 12 GeV(0.1/∆) for MATHUSLA, while p geo = m 2 A − 4m 2 χ /2θ with a typical angle θ where the detector is located: θ 0.5 for MATHUSLA and θ 2 × 10 −3 for FASER. Meanwhile, we take p max to be the same as what we have taken in the text: p max = 100 GeV for MATHUSLA and p max = 7 TeV for FASER.
To determine the reference values in eqs. (4.7) and (4.9), we take m χ = 2.0 GeV, and we read the boundary values of at m χ = 2.0 GeV from the sensitivity plots with ∆ = 0.1, which are shown as green and cyan dots in the figure. In the case of MATHUSLA, we obtain the reference values defined in eqs. (4.7) and (4.9) from this calibration as follows. Here, (cτ χ p max m −1 χ ) 0 , (cτ χ p min m −1 χ ) 0 are the reference values d 0 (p max ) and d 0 (p min ) in the formulae, respectively. On the other hand, in the case of FASER, we obtain the reference values from this calibration as follows.
We show the boundaries with our approximation formulae as black lines in the figure. In both panels, their line types are the same as the sensitivity plots: ∆ = 0.1 (solid), ∆ = 0.05 (dashed), and ∆ = 0.03 (dotted). We just use the reference points on the sensitivity plot for ∆ = 0.1, the other shape of the sensitivity plots is almost reproduced by our formulae even for other choice of ∆. Since we assume that χ 2 is produced only via the Drell-Yan process, the fitting of lower sensitivity curves gets worse for m χ 1 GeV, where the Bremsstrahlung process would be more important.

C.2 SeaQuest
We show the validity of our fitting for the SeaQuest sensitivities on the decay of dark vector mesons in figure 8 by applying our approximation formulae, eqs.  [164] for SeaQuest: for fixed m π /f π = 3 (left) m π /f π = 4π (right). The magenta lines show their results, and black and red lines show our fitting: magenta thick lines correspond to the sensitivity lines for visible decay of dark vector meson, 3-body decay (V →π + − ) and 2-body decay (V → + − ), while magenta thin lines correspond to the sensitivity lines only for 2-body decay. We take the dark pion mass to be m π = m A /3 and the mass of the dark vector meson to be m V = 0.6m A . We take the reference points (magenta dots in the right panel) for our fitting.
ref. [164]. In the literature, the dark photon is heavier than the dark hadrons. The on-shell A is produced through the SM meson decay, Bremsstrahlung, and Drell-Yan at SeaQuest. The produced dark photon mainly decays into the dark mesons, in particular A →ππ and A →πV where V andπ collectively denotes vector mesons and pseudo-scalar mesons, respectively.
The U(1) D -neutral dark vector mesons mix with dark photons as with the ρ 0 -γ mixing in the SM, and then the neutral mesons decay into the SM lepton pair, + − , through the dark photon mixing. On the other hand, the U(1) D -charged dark vector mesons decay into the SM leptons through the off-shell dark photon by emitting U(1) D -charged dark pions. The decay length of the three-body decay tends to be longer than that of the two-body decay due to three-body phase space factor.
Magenta lines in figure 8 show sensitivities of the dark vector meson V decay at SeaQuest: thick lines correspond to 3-body decay (V →π + − ) and 2-body decay (V → + − ), while thin lines correspond to the sensitivity lines only for 2-body decay. The thick sensitivity lines are composed of two parts: upper parts from the three-body decay Γ(V → π + ), while lower parts from the two-body decay Γ(V → + ). It is assumed the ratio of the dark pion mass and the pion decay constant to be m π /f π = 3 (= 4π) in the left (right) panel. In both panels, the dark pion mass is m π = m A /3 and the mass of the dark vector meson is m V = 0.6m A . The produced numbers of V depend on m π /f π , and increases as m π /f π gets larger. In particular, the branching fraction of V is almost unity when we take m π /f π = 4π. The reference values in eqs. (4.18) and (4.19) are determined by the point at m A = 0.1 GeV (magenta dots in the right panel). We use the plot for m π /f π = 4π for the calibration, and hence the branching fraction to the vector meson is assumed to be unity. The lower boundary of the three-body decay is hidden since the boundary is overlapped with the sensitivity to the two-body decay. Therefore, we take the reference points from the upper boundary of the three-body decay and the lower boundary of the two-body decay. We obtain the reference values as follows. Here, [cτ V (V →π + − )p max m −1 V ] 0 and [cτ V (V → + − )p min m −1 V ] 0 correspond to the boosted decay lengths referred as to d 0 (p max ) and d 0 (p min ) in the text, respectively. The black lines show the boundaries computed by our approximation formulae, and the boundaries well coincide with the sensitivity computed by ref. [164]. Since we do not include the dark photon production through Drell-Yan, the boundaries are cut at a few GeV. The black-dashed lines in the left panel also show the boundaries with our formulae. Although the upper boundary of the sensitivity is not shown below 1 GeV in the literature, it is expected that our formulae well explain its extrapolation. By use of our formulae, we also evaluate the upper-boundary of the two-body decay and lower boundary of the three-body decay, which as shown as red-dashed lines in figure 8. The sensitivity of the upper boundaries is also computed by ref. [164], which is shown as the thin-magenta lines, and our formula roughly reproduce their result.

C.3 E137
As for the existing constraint on the dark pion decay, we use the constraint on the visible dark photon decay by E137. We take the reference point at m A = 0.1 GeV, and read out the kinetic mixing on the upper and lower boundaries of the sensitivity area. Figure 9 shows our fitting of dark photon visible decay with the constraint by ref. [4]. The black lines show our fitting, dashed light green line shows the boundary of the existing constraint by ref. [4], and the green points are our reference points. At these points, we get the reference values as follows. Our fitting works well for the dark photon visible decay search at E137.

D Production and geometric efficiency
We use the produced number of the long-lived particles given by eqs. (4.5) and (4.14) for the visible decay search at the LHC lifetime frontier. Here, we naïvely assume that dark nucleons and dark pions are produced when off-shell dark photons with a fictitious mass of the dark dynamical scale are produced. This produced number depends on multiplicity of the produced particles, and we assume n p = 0.04 and n π = 2.0 in the text, which are the analogy of the SM particle production at the J/ψ threshold [155]. Since pions are expected to be more produced than nucleons when they are much lighter than nucleons, the multiplicity of pions is larger than that of nucleons In our study, the mass difference of dark nucleons and dark pions can be smaller than that of the SM hadrons. Therefore, we change the multiplicity of hadrons: in particular, to be the similar multiplicity, n N 1 = 0.1 and n π = 0.2. We show sensitivity plots in figure 10 when the dark sector hadrons are produced on an almost equal footing. In this figure, we take sin θ V = 5 × 10 −3 and ∆ N = 0.1 for the visible decay of dark nucleons, while m π /m A = 1.9 for the visible decay of dark pions. We commonly use U(1) D coupling to be α = 0.05 (α = 7 × 10 −3 ) for n g = 1 (n g = 8). Figure 10 show the changes of the sensitivity plots of FASER (left panel) and MATHUSLA (right panel), and both panels include the visible decays from dark nucleons and dark pions. The color codes are the same as the previous figures for FASER (MATHUSLA): the blue-solid (green-solid) lines for m N 2 = 8.5 GeV and the cyan-solid (yellow-solid) lines for f π = 0.8 GeV, while the blue-dashed (green-dashed) lines for m N 2 = 1.1 GeV and the cyan-dashed (yellow-dashed) lines for f π = 0.1 GeV. Thick lines show the sensitivity curves with eqs. (4.5) and (4.14), while thin lines show the sensitivity curves with different choice of multiplicities. The pion multiplicity smaller than unity implies only three dark pions are required as the final state through the dark hadronization at one hadronic event. Hence, the dark-pion sensitivity curves extends to 3m π < m N 2 in the case of n π = 0.2.
We show the parameter dependences of dark nucleon sensitivities at the LHC lifetime frontier in figure 11. Existing constraints by dark photon searches and the direct detection are the same as in figure 4. For the larger ∆ N , the decay with on-shell dark photon easily -57 - For the latter case, the dark pion sensitivities extend to 3m π < m N2 since the multiplicity is below unity. As in figures 4 and 5, we also show the existence constraints and the future sensitivities to the visible dark photon decay (the top and left-bottom shaded areas and the future sensitivity curves from Belle-II, LHCb, SHiP, and SeaQuest).

JHEP03(2022)176
opens, and therefore the sensitivity to the three-body decay of the nucleons gets narrow in m A axis. However, it is sensitive to the smaller since the lifetime gets shorter. Due to the same reason as what we explained in the text, the sensitivity range gets shrink once we take tha larger couplings or the larger DM mass. Figure 11 shows that some sensitivity curves disappear when sin θ V = 5 × 10 −2 and ∆ N = 0.3.
As for the dark pions at SeaQuest, they are produced via a virtual dark photon with the fictitious mass of order of the SM ρ mesons, and we use eq. (4.17) to estimate their produced number. The proton Bremsstrahlung dominates their production, and the main contribution comes from the ρ meson resonances due to the broad width of ρ mesons. As for the dark photon searches at the SeaQuest [21], as we discussed in the text, the dark photons from meson decays have less energetic and are produced with a larger angle off from beam-axis in comparison with the dark photon from proton bremsstrahlung. Therefore, the geometric efficiency at SeaQuest decreases for the dark photons from meson decay. Similarly to the dark photons, the searches for dark pions can be affected by the geometric efficiency. Even though the production of dark pions is basically from proton Bremsstrahlung, the geometric efficiency for dark pions could change since ρ meson dominates their production. In other words, the invariant mass of final states, √ s f , is close to the ρ meson mass, √ s f ∼ m ρ , and hence the efficiency could be reduced in analogy to the efficiency of dark photons from meson decay. Figure 12 shows the sensitivities of SeaQuest to the dark pion searches, including variation of the geometric efficiency at SeaQuest. The different line-types correspond to the different values of m π /m A : m π /m A = 1.9 (solid) and m π /m A = 1.3 (dashed). We . We take f π = 0.8 GeV (left) and f π = 0.1 GeV (middle) for the sensitivity of SeaQuest, while f π = 0.1 GeV for the exclusion limit from E137 (right). We take the different values of m π /m A : m π /m A = 1.9 (solid boundaries) and m π /m A = 1.3 (dashed boundaries). For the SeaQuest plots, we take the production part to be ten times smaller than that used in the text, which is illustrated as thin lines in the figure. For the E137 plot, we again take the production part to be ten times smaller than that used in the text, which is illustrated as dashed-black lines in the figure.
take different decay constants, f π = 0.8 GeV (left) and f π = 0.1 GeV (right). When the geometric efficiency is ten times smaller than the Bremsstrahlung, the sensitivity areas get smaller and the sensitivity lines are shown as thin lines in the figures.
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.